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ABSTRACT 

The Lyman-alpha forest has become an important tool for measuring the mass power spectrum at high 
redshifts (z = 2 — 4). A crucial intermediate step is the measurement of the transmission power spectrum. 
We present new methods to minimize the systematic and random errors for such a measurement, and 
discuss their implications for observing strategies. Sources of systematic errors explored include metal 
line contamination and continuum-fitting. We advocate the technique of trend-removal in place of 
traditional continuum-fitting - here, a spectrum is normalized by its (smoothly varying) mean rather 
than its continuum - this method is easily automated and removes biases introduced by continuum- 
fitting. Moreover, trend-removal can be easily applied to spectra where continuum-fitting is difficult, 
such as when the resolution or signal-to-noise is low, or for spectra at high redshifts. We further show 
that a measurement of the continuum power spectrum (plus a related quantity) using trend-removal, 
from either low redshift quasar spectra or the red-side of Lyman-alpha, can be used to constrain the 
amount of spurious large scale power introduced by the uncertain continuum, and in principle allows 
the removal of such contamination and thereby expanding scales probed to larger ones. We also derive 
expressions for the shot-noise bias and variance of the power spectrum estimate, taking into account 
the non-Poissonian nature of the shot-noise and the non-Gaussianity of the cosmic fluctuations. An 
appropriate minimum-variance weighting of the data is given. Finally, we give practical suggestions on 
observing strategy: the desired resolution and S/N for different purposes and instruments, and how to 
distribute one's finite observing time among quasar targets. Also discussed is the quasar spectroscopic 
study of the Sloan Digital Sky Survey (SDSS), which has the potential to measure the power spectrum 
accurate to better than 1% per mode (Ak ~ lCP 4 s/km) - the techniques presented here will be useful 
for tackling the anticipated issues of shot-noise and continuum contamination. 

1. INTRODUCTION 

Recent theoretical research on the low column density (A%i < 10 145 cm~ 2 ) Lyman-alpha (Lya) forest at redshifts z ~ 2—4 
points towar ds a picture in which the forest consists largely of mildly nonlinear fluctuations of a smooth intcrgalactic 
medium (e.g. Bi ct al. 1992|; Gen et al 1994 ; Zhang ct al. 1995; iReis ene; 
Escuclc et al. 1996j; |Mucckct et al. 1996|; |Bi fc Davidsen 



Miralda- 



1997 



er fc Mira.lda-Escudc 1995 1 Hemquist ct al. 1996 



Bond fc Wadsley 1997| ; |Croft ct al. 1 
This provides the motivation to 



997|; Hguj 

) analyze 



et al. W\) t[ |tiui &: Uncdm liJiJY| ; see |Kaucn lor a review and lurtner ret. J. Tins provides trie motivation to analyze 
the quasar (QSO) absorption spectrum as a continuous field with fluctuations, rather than as a collection of discrete 
absorption lines. The two-point correlation or its fourier transform, the power spectrum, comes to mind as a useful 
and common statistic used in other areas such as the microwave backgro und or galaxy large-scale-structure. Indeed, its 



application to QSO spectrum has bee n discussed by a number of authors (Zuo fc Bond 1994; Miralda-Escude et al. 1996| 



Bi fc Davidsen 1997| ; |Cen et al. 1998[ ). Recently, Croft et al. ( |1998a|; |1998bp (sec also |Hui 1999] ; |McDonald et al. 1999j 
have shown that the mass power spectrum can be recovered fr om the QSO transmi s sion power spectrum , from which 
one could further deduce cosmological parameters such as fl m (|Wcinberg ct a l . 199S ; Phillips et al. 200C ). There exist 



at present a large number of high quality QSO spectr a (e.g. Hu et al. 1995| : |Lu ct al. 1996 ; Kirkman & Tytler 1997 



Cristiani et al. 1997; Kim et al. 1997; Rauch et al. 1997) which makes this an exciting field of research. Upcoming quasar 
surveys such as the Sloan Digital Sky Survey (SDSS) and the Anglo- Australian Telescope Two Degree Field (AAT2DF) 
will enlarge the database significantly. 

Here we take the view that the QSO transmission power spectrum / correlation is interesting in its own right, and focus 
on how to best measure it from the observed QSO spectra, independent of the underlying physical picture of the forest. 
The two major questions are: (1) what are the main sources of systematic errors and what are the best ways to bring 

1 Institute for Advanced Study, Olden Lane, Princeton, NJ 08540; NASA/Fermilab Astrophysics Center, Fermi National Accelerator Labo- 
ratory, Batavia, IL 60510; e-mail: lhui@ias.edu 

2 Dcpartmcnt of Astronomy and Astrophysics, University of Chicago, Chicago, IL 60637 
3 Physics Department, Princeton University, Princeton, NJ 08544 
4 Department of Astronomy, California Institute of Technology, Pasedena, CA 91125 
5 Canada- France-Hawaii Telescope, 65-1238, Mamalahoa Hwy, Kamuela, HI 96743 

6 Center for Astrophysics and Space Sciences, University of California, San Diego, MS 0424, La Jolla, CA 92093-0424 

1 



2 



them under control? (2) how to estimate the shot-noise, and to best combine data with different signal-to-noise (random 
errors)? 

Let us start by defining the transmission power spectrum and correlation function. Two possibilities arise. One of 
them we call the un-n ormalized power spectrum P un / two-point correlation £ un (Weinberg f 998, private communication; 



McDonald et al. 199S) 



= (f(u')f(u' + u)) , P un (k) = J U{u)e- lku du (f ) 

where / is the transmission defined by / = e~ r with r being the optical depth (the absorption is then 1 — /), Horn' 
is the observed velocity (or redshift or wavelength) along a line of sight, and k is its fourier counterpart. The angular 
brackets denote ensemble averaging. 



The other we call the normalized power spectrum P / two-point correlation £ (e.g. Zuo & Bond 1994) 



f(u) = (S f (u')5 f (u' + u)) , P(k) = / £(u)e- ikn du (2) 



where Sf = (/ — /)//, with / being the mean transmission (/}. 

We will almost exclusively focus on the latter, but will discuss at some point the pros and cons of the two, especially 
with regard to systematic errors. Unless otherwise stated, hereafter, power spectrum / correlation refers to the normalized 
version. 

The layout of the paper is as follows. In §0, we provide a brief overview of how the raw data output (a two-dimensional 
CCD image) is reduced to a one-dimensional QSO spectrum. Note that the quantity / above is never observed directly. 
It is important to have a description of how the whole data reduction procedure works, which is sometimes hard to find 
in the literature. We give some illustrations by showing simulated spectra with various realistic levels of noise. 

In §Q, we discuss the esti mat ion of the power spectrum and two-point correlation, beginning with the introduction of 



the quadratic estimator in §4.1. An important point pertaining to the estimation of the two-point correlation is raised 
here - most estimators employ ed in the literature ar e not optimal; an alternative is given here which is an analogue 



of the Landy-Szalay estimator (Landy & Szalay 1993) introduced or igin ally for galaxy surveys. Aside from this point, 
we focus exclusively on the estimation of the power spectrum. In §]4.2| , we discuss three sources of systematic errors: 
continuum-fi tting, gaps and metal absorption lines. Particular attention is paid to issues related to continuum- fitting. We 
advocate in £|4.2.l| trend-removal to replace traditional continuum-fitting, which avoids the latter's pitfalls. We further 



propose in § [4.2.2| that the power spectrum of the continuum can be estimated using trend-removal as well, which offers 
us a way to measu re acc urately the transmission power spectrum on large scales where the continuum fluctuations mig ht 



be important. In §4.2.3, we discuss the effects of gaps and un-removed metal lines. We then turn our attention in §4.3 to 
random errors. We emphasize here that the shot noise is not exactly Poisson distributed because of the particular way the 
data is reduced. We point out the importance of subtracting the shot-noise-bias correctly, and describe a systematic way 
of assigning error-bars to the power spectrum, and introduce minimum-variance-wcighting techniques to combine data 
of different qualities. Some results here are stated without justification. The aim in this section is to summarize useful 
results for readers who might not be interested in details of the derivations, which are provided in the Appendices. The 
techniques used in the Appendices should be of broad interest e.g. the issue of generally non-Poissonian shot-noise might 
be relevant for galaxy power spectrum estimation. 

Lastly, we conclude in §^| We summarize here our recipe for transmission power estimation - readers who would like 
a quick overview of our methods can skip directly to this section, and only refer back to the relevant sections to fill in 
the details. We give general advice on observing strategies, and discuss in particular analysis issues relevant for the Sloan 
Digital Sky Survey. 

Before we begin, let us first make a few clarifying remarks about some of our notation and terminology. 

2. TERMINOLOGY AND NOTATION: AVERAGING AND AVERAGES 

In this paper, we refer to two different kinds of fluctuations which should be clearly distinguished. Take the observed 
photon (or electron) count from a quasar as an example, Nq. As one moves along a spectrum, the photon count fluctuates 
because of two very different reasons. First, it fluctuates because the universe is intrinsically inhomogeneous, giving rise to 
non-uniform absorption - we will refer to these as cosmic fluctuations. Second, it fluctuates because the observed photon 
count is a discrete realization of the underlying cosmic signal - Poisson fluctuation is the canonical example but not the 
only possible one, we will refer to these as discrete fluctuations. 

We define two different kind of averages corresponding to these two different kind of fluctuations. The discrete average 
of the observed photon count is denoted by (Nq)d = Nq. In other words, Nq constitutes a discrete realization of the 
underlying cosmic signal Nq. This signal Nq itself suffers from cosmic fluctuations, and we will denote its ensemble 
average by (Nq) — Nq. A fixed quasar continuum is assumed in this ensemble average i.e. it is the fluctuation in the 
spectrum caused by intervening absorption that constitutes the cosmic signal we are after. Finally, we will sometimes use 
() to implicitly stand for (()d) e.g. (Nq) actually means ((Nq)d) which is the same as (Nq) = Nq [] 

To recap: 



7 Two exceptions in the use of (): in § |4,2.2 , we use () to include averaging over an ensemble of different continua, and in § [l.3| we use ()kk 



to denote averaging over a shell of Fourier modes 
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• Nq is the directly observed quasar photon count. 

• Nq = (Nq)d refers to the idealized quasar photon count if one has infinite signal-to-noise (S/N). 

• Nq = (Nq) = ((Nq)d) refers to the quasar photon count if one has infinite S/N and if one averages over all possible 
cosmic fluctuations keeping the continuum fixed e.g. by taking the same quasar and putting it at all possible 
orientations in the sky. For instance, if ./Vq = Nce~ T where Nq is the true continuum and r is the optical depth, 
then Nq = Nc(e~ T ) where (e~ T ) is the mean transmission. 

Note that, when we use the term discrete, it is not implied that Nq is necessarily an integer, although it is derived 
from some integral quantity such as the electron/photon count. We use the term discrete fluctuations instead of the usual 
Poisson fluctuations, because as we will see, Nq is often not Poisson-distributed i.e. (Nq)d — (Nq)^ ^ ^Vq = (Nq)d 



(see §3.1). The term shot-noise is often used to describe Poissionian discrete fluctuations, but in this paper we will use it 
more broadly to include non-Poissonian discrete fluctuations as well. 

Finally, we note that we use the term 'quasar counts' throughout this paper to refer to the photon counts from a quasar, 
rather than the number of quasars in a given patch of sky. 

3. DATA REDUCTION: FROM THE CCD IMAGE TO THE QSO SPECTRUM 
3.1. A Brief Description 

We discuss briefly here aspects of the data-proces sing n ecessary for un dersta nding the noise prop erties of the reduced 
quasa r spectrum. The read er is referred to Home ( 1986 ), Zuo & Bond ( 1994 ), Barlow & Sargent ( 1997| ), Rauch et al. 



(1997) and Cen et al. (1998) for more discussions. 

The raw data consist of a two-dimensional (spatial and spectral) array of counts (data numbers or photon counts 
converted from them) from a CCD image. The one-dimensional array of estimated quasar counts in the spectral direction 
(as a function of velocity, redshift or wavelength) is obtained by collapsing the data in the spatial direction in the following 
fashion: 

^ = ^rV(^ w -^) (3) 

We have introduced and will stick with the following notations: the Latin index such as i and the Greek index such 
as a denote the spatial and spectral coordinates respectively of a CCD pixel (there are in fact a few exceptions, which 
should be clear from the context); Nq is our estimated quasar count, N^ AW is the raw count, and N^ is the mean 

background count which includes the sky and the readout offset; W l13 is a weighting of the spatial pixels for each spectral 
coordinate (3, and W a ^ represents a rebinning of the spectral pixels that is sometimes done to achieve, for instance, a 
linear wavelength scale. Note that the a and i dimensions do not necessarily align with the two perpendicular axes of the 
CCD chip. The optical setup could be such that the spectrograph slit appears tilted at an angle to the CCD axes. We 
use " to emphasize the fact that the quantity of interest is a random variable with fluctuations. The ~ denotes a discrete 
average: e.g. N^ = (N^)d where ( )rj denotes discrete averaging. 

Implicitly assumed in the above formulation is that the discrete average N^ is known, which is of course not strictly 
true, but since a typical slit covers a significant number of pixels that do not have any quasar photons, and since the 
background is often quite uniform, the discrete average can be estimated to high accuracy. Note also the readout offset 
can be measured using short exposures with closed shutters or from the CCD overscan region. 

The weighting W 1 ^ typically has non-trivial spectral dependence (J3) to remove at least two artifacts: variations in 
detector efficiency across the chip and a non-flat blaze. The former is usually estimated in a procedure called flat-fielding 
by shining a lamp into the detection system. The latter arises because of the non-trivial shape of a diffraction order. 
This can be partially estimated in the flat-fielding procedure, but is best done using a spectrophotometric standard star, 
usually a white dwarf. While the correction for the first artifact should be quite accurate, the blaze-removal is often 
approximate. Any residual that is not correctly removed will show up in the form of a non-trivial effective continuum. 
We will see in §|4.2| perhaps some evidence of it. We assume in this paper that such artifacts show up as fluctuations on 
large scales (since the blaze itself is smoothly and slowly varying across a given order) but not on small scales (we will 
quantify the scales later on). 

To make the above concrete, the raw count is given by: 

<w=^ + < (4) 
where the quasar contribution has the following discrete average: 

NS = {NS)D=g^N^ (5) 



where g l ^ B is the point-spread function which describes how the light from the quasar gets spread-out in the spatial direction 



i at a given spectral coordinate f3, and g£ accounts for the variation of the blaze and quantum efficiency as a function of 



wavelength. The symbol Nq denotes the underlying quasar count (or cosmic signal i.e. discrete averaged). 
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Many different rebinning kernels W a @ (eq. |3)) are possible. The simplest choice is of course no rebinning with 

W a p = 5 af3 . 

There are several possible choices of the weig hting W lf3 (eq. H), but any sensible choice has to satisfy the requirement 
that (^Vq)d = Nq, up to some constant normalization factor. This assumes that artificial fluctuations introduced by the 
blaze or detector efficiency are correctly taken out. If not, it shows up effectively as part of the continuum. 

We give two examples of W l P here. The first is basically a uniform weighting over the spatial pixels that correspond to 
a given spectral coordinate: 

W^ = 1/G£E<#)' (6) 

i 

where the range of i, or j, is chosen to lie within, say, some fraction of the quasar seeing disk. There is sometimes an 
additional complication due to cosmic ray hits, which will be discussed below. 

The seco nd is a minimum variance weighting (different fr om m inimum variance weighting for measuring the power 
spectrum; §4.3) over the spatial pixels, introduced by Home ( [19861) : 



W i0 = (l/^)(^/Kiw)/(E[3ps] 2 /^lw) (7) 



where is the variance in the raw count: 



ViU = ((N^ w f) D - (N^ w f D = ivf + Vf (8) 
V B sf> = iVf + v£? . 



where Vg , the background variance, has two contributions, the sky variance iVg and the readout variance V^ q . A 

rjp 

'RAW 



word of caution is necessary here regarding the second weighting. The raw variance, depends on the underlying 



cosmic signal or quasar count {Nq i.e. discrete averaged) which is not directly observable (the discrete averaged sky 
count and the true readout variance are also strictly speaking not directly observable, but they can be estimated quite 
accurately because they are relatively uniform and can be observed over a larger number of pixels). Modeling Vp^w 
using the measured raw count (i.e. using N Q P instead of N^ in eq. (|]) could lead to a biased estimation of Nq. Home 



(1986) suggested an iterative scheme to avoid this problem, but implementations of this weighting should be checked for 
a possible bias. 

Pixels affected by cosmic-ray hits, which are usually easy to identify because of their wild fluctuations and spiky nature, 
are dealt with in two different ways, depending on the severity. For a given spectral coordinate, if only a small fraction 
of the corresponding spatial pixels are affected, the weighting in eq. (g) or eq. (R) is simply modified by allowing i and 
j to only range over the unaffected spatial pixels. However, if all or most corresponding spatial pixels are affected, then 
all recorded counts at that spectral coordinate are discarded, leaving a gap in the reduced quasar spectrum. Gaps could 
result also because of metal-line removal (an alternative would be to fit for the metal- line and subtract, instead of simply 
discarding the pixels) or defects in the CCD. 

Finally, the (random) error-array output at the end of the data reduction corresponds to an estimate of 



((N%-N$)*)d = /E^^Xiw (9) 

V i,P 



the error-array is only an estimate of the above quantity, because the true V^f AW is unknown, but is estimated using the 



where we have assumed that the noise-fluctuations are independent among the pixels. We emphasize that in practice 
the error-array is only an estimate of the above quantity, t 
observed raw counts (using Nq instead of Nq in eq. ||). 

It is clear from the above discussion that in general fluctuations in Nq are non-Poissonian, in the sense that (Nq)d — 
(Nq) 1 ^ ^ (Nq)d — Nq. This is because of several reasons. First, Nq suffers from additional discrete fluctuations from 
the background counts. Second, the weig hting W a(3 and W i0 are in general non-trivial (i.e. not unity). A very simple 
example of the effect of non-unit weights is: suppose we multiply a Poisson variable by a factor of 2 and call the result y, 
it is easy to see that (y 2 ) D - (y) 2 D = 2(y) D ^ (y) D . 

For the rest of this paper, we will pick for simplicity the weighting kernels W a /3 — 6 a /3 and Wip as given by eq. (@) . In 
the Appendices we will indicate where some of our expressions have to be modified to account for more general weightings. 

3.2. Simulated QSO Spectra 

For illustrations, and for later analyses, we have generated several different QSO spectra. The underlying noiseless 
(theoretical) transmission (/ = e~ r ) is shown in the bottom panel of Fig. |Tj, and its associated power spectrum is shown 
in the to p pan el of the same figure. This is drawn from a SCDM (Standard Cold Dark Matte r) simulation discussed in 
Gnedin ( 199S| ) which made use of the Hydro-PM algorithm developed by Gnedin & Hui ( |1998| ). The cell-size (comoving 
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size of 10h _1 kpc) is small enough to resolve the effective Jeans scale, and so should retain all small scale structures. 
However, the box-size is unrealistically small (comoving size 2.56 h _1 Mpc) which means a significant amount of large scale 
power is missing. For most of our investigations here, it is not necessary that the simulations are highly realistic, but our 
simulated transmission power spectrum is in fact broadly consistent with an observed one (Fig. HJ). The long line of 
sight in Fig. [l] is generated by shooting a ray at some oblique angle through the simulation box and allowing it to wrap 
around the box several times , but never repea ting itself. The mean redshift here is z = 2.85. The ionizing background is 



chosen to give (e T ) = 0.64 ( Press ct al. 1993 ) 



An example of a somewhat realistic reduced QSO spectrum (Nq in eq. |3J) and its error array (eq. |9J) can be found 
in the bottom two pan els of Fig. ^ (ignore the other two panels for the moment). They are generated based on the 
prescriptions given in §3.1, assuming W a @ — 8 a ^ and W 1 ^ is given by eq. (||). Briefly speaking, what we do is to 



first generate an array of g which represents g£ ^ . g^" (i.e. we do not actually simulate the full two-dimensional CCD 
image; the spatial dimension is collapsed into g a ); then, we create a Poisson realization of the (intermediate) quasar count 
g a N^,e~ Ta where Nq is the continuum and e~ Ta is predicted by our cosmological model; we similarly create a Poisson 
realization of the background count g a x const, where the constant represents some fraction less than 1 and then subtract 
from it its Poisson mean, the end-result is then added to the above quasar count; lastly, we divide by g a to obtain the 
reduced quasar count Nq. 

Note that the overall level of the reduced quasar count can be scaled up or down (because we are not interested in 
the absolute brightness of the quasar), provided the error array is scaled accordingly to conserve signal-to-noise (S/N; 
count divided by square root of the variance) . This example resembles a high quality Keck HIRES spectrum, with S /N 
reaching up to 100 at certain pixels. It is composed of 12 echelle orders, 50A each (e.g. the instrument HIRES on the Keck 
telesc ope is an echelle spectrograph which consists of two diffraction gratings crossed at 90° to each other; see Vogt et al. 



1994|) . The pixel size is 0.05,4 with a resolution Full- Width-Half-Maximum (FWHM) of 0.125A The example represents 
a case in which a relative calibration (but not necessarily absolute fluxing) between the orders has been attempted. The 
dashed line in the bottom panel shows the input continuum. 

The error-array shows a lot of variations. A model of the random error as Gaussian distributed with uniform S/N that 
is sometimes used in the literature misses much structure. About 3% of the spectrum consists of gaps which arise due 
to severe cosmic-ray hits. The spikes in the error-array correspond to wavelengths at which a fraction, but not all, of 
the corresponding spatial pixels are affected by cosmic-ray hits. They also take up 3% of the spectrum. It is easy to see 
how these spikes arise from eq. @ and @. At wavelengths where some of the spatial pixels are thrown out because of 
cosmic-ray hits, Wip is enhanced because the sum over j in its denominator is restricted to fewer pixels. Since it is the 
square of Wip that enters into the variance, a modest enhancement becomes a spike. Note how for each echelle order, 
the S/N drops towards the two ends. This is because of the blaze function which tends to suppress the flux at the ends. 
Note also that the S/N has a general decline towards the blue. This is due to a combination of a falling continuum, and 
decreasing detector efficiency. 

Sometimes, a relative calibration between echelle orders is either difficult or simply not attempted. An example is 
shown in Fig. [| Note how the continuum is broken into 12 discontinuous pieces]^] These are taken from continuum-fits 
to actual data. 

An example of data with much poorer quality is shown in Fig. [?]. The pixel size is 0.5A and the FWHM is 1.17A. The 
S/N is about 10 times worse than the two examples above. Such a spectrum could be the output of, say, a low dispersion 
single-grating spectrograph, which does not have the characteristic division into short pieces as in the case of the echelle 
spectrograph. 

All other simulated data in this paper are slight variations of the above, which will be described in turn at the appropriate 
places. 

4. ESTIMATING THE POWER SPECTRUM / TWO-POINT CORRELATION 
4.1. The Quadratic Estimator 

Given the one-dimensional array of estimated quasar counts Nq (eq. how should one go about estimating the 

two-point correlation or the power spectrum? 

A common practice is to first continuum-fit, i.e. to estimate the continuum level Nq and divide JVq by it to obtain an 
estimate of the transmission / = Nq/Nq. Then, the estimators for the un-normalized two-point correlation and power 
spectrum (eq. [0) are: 



a, (3 a,0 



where b(k) subtracts out the shot-noise (i.e. a bias), and w a ^(u) and w a ^(k) are weighting kernels for which we will give 
some examples shortly (to be distinguished from W a @ in eq. H). These are commonly called quadratic estimators for 
the simple fact that they make use of the_data f a in quadratic combinations. 

Defining the mean transmission to be /, the obvious extensions of the above estimators, for the normalized two-point 
correlation and power spectrum (eq. Q; unless otherwise stated, the two-point correlation or power spectrum with no 

8 In some cases where the different echelle orders overlap, there could be two jumps at each order junction. 
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qualifications refer to the normalized version), are 



£ w^{u){t - f){P - DIP , Pi(k) = J2 ^(k)(f a ~ f)(f P ~ 1)1 P ~ b(k) 



(11) 



a, (3 



a, 



However, the form of the power spectrum or two-point correlation estimator given above suggests an interesting variation 
which allows us to avoid continuum-fitting altogether: (f a — /)// can be estimated instead by (N a — N a )/N a where 
N a is the mean count defined by N a = (N a ). Here ( ) denotes the cosmic average i.e. this corresponds to averaging out 
the cosmic fluctuations in /, for a fixed continuum ((N a ) = Ngf where Ng is the true continuum). [] Note that the 
mean count is dependent upon a because of the slowly varying continuum. We will discuss how to estimate N a shortly. 
The key here is that the absolute level of the continuum gets divided out by definition. Hence, let us define the following 
alternative estimators of the two-point correlation and power spectrum: 



a, (3 a, (3 

51 = (N% - N$)/N% 



(12a) 



(12b) 



This alternative power spectrum estimator is what we will focus on, but we will also briefly investigate the behavior of 
the estimators in eq. (|lo| ) and ([ll|). 

It remains to specify what w^(u), w a @(k) and b(k) are. The simplest choice is to use uniform weighting i.e. for the 
two-point correlation, it corresponds to: 



w 



a/3 



(13) 



where 9 a ^(«) is equal to one if the two pixels a and (3 are separated by a distance u (or more generally, the distance falls 
into a bin that is centered around u with some finite width), and zero otherwise. Using the above w a ^(u) corresponds to 
simply counting all pairs separated by a distance u, normalized by the total number of pairs. 

With the above weighting, eq. ( |l2|) is analogous to an estimator of the two-point correlation introduced by Landy 
& Szalay (|l993|) for galaxy surveys: (DD - 2DR + RR)/RR, if one identifies DD with ^^^(iijjV"^, DR with 



w a P(u)N a NP and RR with J £2 a pW af3 (u)N c "N 13 and assumes N a varies very slowly with a on the scale of interest 



u (th e analogy becomes exact in the limit of a uniform N a ). As shown by Landy & Szalay ( 1993| ) (see also [gzapudi "fc 
Szalaj 199^; Podelson ct al. 1997D , a common alternative estimator DD/ RR — 1 (equivalent to the estimator used by e.g. 
Zuo & Bond 1994: Ccn ct al. 199§| ) is actually less desirable as it gives a larger variance compared to (DD — 2DR+RR)/ RR. 

With this being said, we are going to focus our attention on the power spectrum from now on, although most of our 
treatments below can be applied to the two-point correlation as well. For the power spectrum, the simplest choice of 
uniform weighting corresponds to: 



w ap {k) = {C/M 2 )R aP (k) , R at \k) = (l/n fe )^e- lfe > Q -"^ 



b(k) 



E 



VS 



( n q) 2 



(14a) 

(14b) 



where C is the total length of the spectrum (in whatever units one prefers) and TV is the total number of spectral pixels, 
R a °(k) represents an average of the Fourier basis over some bin or band in k-space i.e. we estimate the power spectrum 
at a bin centered at k by averaging over contributions from each k' that belongs to the bin (jik is the total number of 
modes in it). ^ The symbol b(k) represents the shot-noise contribution to the power that has to be subtracted off, 
and represent the point-spread function and the blaze as in eq. (||), and V§ is the background contribution to the 
shot-noise (including the sky and readout, see eq. ||). Note how a part of the shot-noise depends on the reciprocal of 
the mean quasar count, as expected for Poisson fluctuations, quite analogous to the shot-noise of galaxy distributions. [^] 



However, the factor q a , which arises from non-trivial weighting of CCD pixels (W ia ; eq. 



signifies that the shot-noise 



is not strictly Poisson distributed. Moreover, there are extra contributions to the shot-noise from the background counts, 
which are generally absent in the case of galaxy surveys. Derivations of the above statements are given in Appendix A. 
The corresponding power spectrum estimator obeys: 



dk' 

(A(fc))= / ^G(MW) , G(k,k') 



uj3 



w al3 (k)e lk ' {ua - uf3) 



(15) 



Sec 



is implicitly assumed that discrete averaging has been carried out before cosmic averaging i.e. we use () interchangeably with {()d)- 



10 See Seljak (1998) and B ond c t al. (1997) for discussions on precautions one should take on binning 
n See e.g. Feldman et al. (Il994|). 
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where G is a window function that resembles , for k 3> l/£, a delta function centered at k = k' with a width of the order 



of 1/C The normalization of w al ^(k) in cq. (14a) ensures that J dk'G(k, k')/2iT = 1. See Appendix A for a derivation. 

It should be emphasized, however, that the above statements are strictly true only if one ignores uncertainties in the 
mean count Nq - i.e. Nq is not known a priori but is instead estimated from the same data from which one is trying to 



measure correlations. We will discuss this further in § 4.2.2] . It suffices to say here that our results in this section remain 



valid as long as one stays away from scales comparable to the entire length of the quasar spectrum. 

4.2. Systematic Errors 
4.2.1. Continuum- fitting versus Trend-removal 



The power spectrum estimator P2 in eq. (|12j), which we are going to focus most of our attention on, requires an estimate 
of mean count Nq. The mean count is not strictly uniform because of a slowly fluctuating continuum i.e. Nq = Ngf 

where Nq is the continuum and / is the (flat) mean transmission. We assume N a has the following form: 

^ a = E c v a (i6) 

a 

where p° is a constant, p 1 is the first order polynomial (p la = u a ), p 2 is the second order polynomial (p 2a = (u a ) 2 ), and so 
on. The coefficients C a need to be estimated from the quasar counts Nq. Note that most of our following arguments would 
go through for a different set of basis functions. One key assumption we will exploit is that N a is slowly fluctuating, so 
that we can truncate the above series at relatively low orders. Continuum-fitting in practice makes the same assumption. 
To estimate C a , we use a linear estimator: 

C a = ^M aa N%, (17) 

a 

where M is a matrix to be specified. Comparing eq. ( ^6|) and Jl7|), it is not hard to see that M has to satisfy J2 a M aa p ba = 
S ab . The simplest choice is to adopt, in vector notation, M = (pp T ) _1 p where pp T in component-form is [pp T ] ab = 
y^ a p aa p ba . In summary, this means our estimator for the mean quasar count is 

N « = Y / La0 K > L = P T (PP T )^ 1 P (18) 



where L in component-form reads L a " = '^2 iah p aa p ah p h ^ with p ah b eing the inverse of the matrix [pp T ] ab . More sophisti- 



cated versions of the above can be found in Rybicki & Press (1992). Our experience is that the simple version given here 
suffices, because the shape of the true quasar continuum is quite uncertain anyway. 

Note the crucial differences between traditional continuum-fitting and an estimation of the mean count as described 
above. The above method makes no reference to the absolute level of the continuum i.e. the count level where there is 
supposedly no absorption. Continuum-fitting in practice often involves human intervention (eye-balling) in the identifi- 
cation of such a level. In contrast, eq. ( fl8h is straightforward to implement and automate. The mean count is then used 
to normalize the quasar count as in eq. ( |l2b| ) before the power spectrum is estimated (eq. [jl2|). We call this procedure 
trend-removal to distinguish it from traditional continuum -fitting. Trend-removal is widely used in other areas (e.g. Press 



et al. 1992a; Rybicki & Press 1992: Tcgmark et al. 1998| ). It is akin to the estimation of, say, the long-term trend of the 



stock market in the midst of all its daily fluctuations. 

Eq. (1^), together with eq. ( |l2| ) and (|14|), completely specifies the main power spectrum estimator we advocate. Several 
tests follow. 
Test 1 

In Fig. ^, we show the effect of different choices of the mean-transmission basis p (eq. |l6|]). The simulated spectrum 
is of Keck- HIRES quality, with a S/N as high as 100 at certain pixels, and it assumes one has a good relative calibration 
between the different echelle orders (12 in all) i.e. an almost ideal, state-of-the-art observed spectrum. The second panel 
from the top shows the recovery of the mean transmission. The solid line represents the true (input) mean. The rest 
shows the recovered mean for different bases: dotted line for a basis consisting of p° only (a constant i.e. the continuum 
or the mean is modeled as completely flat); short-dashed line for a basis consisting of polynomials up to the third order; 
long dashed line for also a basis of polynomials up to the third order but with coefficients fitted separately for each echelle 
order. The short-dashed line seems to give the best match to the true mean. However, none of them is perfect because 
the true mean does not, by choice, have a polynomial shape. This is what is likely to happen in practice - lacking a good 
understanding of the physics that determines the continuum shape of any given quasar, the best one can do is to pick 
a reasonable basis which contains enough freedom to describe the general features of the continuum, but not so much 
freedom that one overfits. 

The important question is what impact the choice of basis has on power spectrum estimation. This is illustrated in the 
top panel of Fig. 0, where the fractional error in the measured power spectrum is shown. The one that gives the best 
match to the true power spectrum is indeed the one where a simple basis of p°, p 3 is used for the whole length of the 
simulated spectrum. The biggest effect of underfitting (dotted line) or overfitting (long dashed) the mean transmission 
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is on the power spectrum estimation on large scales. They cause respectively over- or under-estimation of the large scale 
power spectrum. An additional effect is that overfitting tends to introduce spurious power on small scales as well - witness 
the enhanced fluctuations in the error on small scales for the long-dashed line. We will see this more clearly in the next 
Test. 

Without any prior knowledge of the intrinsic continuum shape of an observed quasar, how does one decide if one is 
overfitting or undcrfitting? One way is to look at the region of the observed spectrum redward of the Lya emission line, 
which is free of the forest, and the continuum is therefore relatively easy to reconstruct. Assuming the general level 
of continuum-fluctuation is the same both redward and blueward of Lya, one can then gain an idea of what a good 
mean-transmission basis might be. Low redshift QSO spectra, where the continuum can be quite easily recovered even 
blueward of Lya, can also be used to gauge the scales at which the continuum fluctuates - naturally, one could also use 
such spectra to check the assumption that continuum-fluctuations have similar characteristics redward and blueward of 
Lya (more on this below). 
Test 2 

In Fig. |^, we took the continuum-fits to an observed quasar spectrum and use them as the input continuum for our 
simulation. The simulated spectrum here represents a case in which no relative calibration between echelle orders have 
been attempted, which is often the case. This is why the continuum in the bottom panel is broken up into 12 pieces. The 
second panel from the top again illustrates the recovery of the mean transmission: dotted line for a flat model-continuum 
for each order, and dashed line for a basis of polynomials up to the third order, also separately for each order. The solid 
line is the true mean transmission. The top panel shows the accuracy of the corresponding power spectrum estimations. 
The assumption of a simple flat continuum for each order gives a power spectrum that is accurate to ~ 1%. On the other 
hand, overfitting with up to third order polynomials not only causes an under-estimation of power on large scales, but 
also creates spurious power on small scales. 

Combining Fig. || & || (note that they show the power spectrum estimation on different scales), the lessons are: 1. it is 
better to err on the side of undcrfitting the mean, which tends to over-estimates the power on large scales, but leave the 
power on small scales relatively unaffected (this relies crucially on the fact that the continuum has fluctuations only on 
large scales); 2. without sufficient prior knowledge of the true shape of the continuum, one can at least make conservative 
statements about the small-scale power, but the large-scale power is likely prone to systematics, unless some correction is 
made. 

One additional comment: the input continuum in Fig. which is taken from fits to actual data, certainly seems to 
sug gest t hat the observed continuum has fluctuations on scales of an echelle order (~ 50 A). (We will quantify this better 
in § |4.2.2 .) ft is unclear whether this is truly due to the intrinsic continuum, or whether it is an artifact of imperfect 
blaze removal or flat-fielding (see §3.1). If it is the former, then ~ 50A represents a fundamental limit beyond which one 
cannot reliably measure the transmission power sp ectrum , at least not without some additional prior knowledge of the 
true continuum (which is what we will discuss in §4.2.2). If it is the latter, then in principle one should be able to do 
better and extend the range of reachable scales to larger ones. Which is the case remains to be seen. 
Test 3 

In Fig. ^, we show the effect of traditional continuum-fitting, which requires some degree of eye-balling. The same 
simulated spectrum as in Fig. |^ is given to an observer (one of the authors) with no knowledge of the input continuum. 
Note that the top second panel now shows the actual continuum level rather than the mean transmission level. The 
estimated continuum actually matches the true one surprisingly well. However, one can still see that the continuum is 
generally underestimated. In the top panel, we show the accuracy of two different power estimates. The long-dashed line 
corresponds to an estimate of the un- normalized power spectrum as defined in eq. (H) (the estimator is eq. |l0[ ) . There is 
clearly a ~ 5% positive bias here because of the underestimation of the continuum. One way to correct for it is of course 
to use simulations: applying exactly the same procedure to the observed data and the simulated data, and see how much 
bias results; but the size of the bias is likely to be model dependent. A simple alternative way to cure this problem is to 
measure the normalized power spectrum instead, using the continuum- fitted data, i.e. using the estimator in eq. (pr|). 
This is shown with a short-dashed line. It has an accuracy of ~ 1%, comparable to the dotted line in Fig. ||. In view of 
this, it seems by-passing continuum-fitting altogether and proceeding simply with trend-removal is desirable. 

Test 4 

The failure of traditional continuum-fitting is more dramatic in cases where there is a lot of absorption e.g. at high 
redshifts. In Fig. |^ is shown a simulated spectrum with the ionizing background adjusted to give a mean transmission of 
0.39, which is about the observed value at z = 4 (Press et al. 1993). The continuum is more seriously under-estimated 
leading to an overestimate of the un-normalized power spectrum by ~ 20% (the upper long-dashed line) . The normalized 
pow er spectrum, estimated either using the continuum-fitted data (eq. |ll| ) or using directly the trend-removed data (eq. 
p2]), is much more accurately measured. 

Test 5 

Another example in which traditional continuum-fitting fails is shown in Fig. [| This is based on the same spectrum as 
in Fig. [|, but convolved with a Gaussian of 1.17 A FWHM and with much poorer S/N compared to the simulated spectra 
above. This is likely not the product of an echelle spectrograph, hence there is no division into 12 orders. We repeat 
the exercise of continuum-fitting and then power spectrum measurement as before. Interestingly, the significant discrete 
fluctuations due to the low S/N here cause an overestimation (unlike in Test 3 and 4) of the continuum level, and so an 
underestimation of the un-normalized power spectrum. Once again, the normalized power spectrum does not suffer from 
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the same problem. Note the somewhat large fluctuations of the estimated power - this is largely due to the high level of 
shot-noise. 

Fig. shows the measurement of power spectrum through trend removal instead. A third order polynomial is used to 
estimate the mean transmission. The resulting (normalized) power spectrum estimate (eq. |l2] ] ) is of comparable accuracy 
to that using the continuum-fitted data. We also show in the top panel as a dotted line the power spectrum estimate 
without shot-noise subtraction (eq. juj). Clearly, shot-noise introduces a bias of the order of 10% here. We will have 
some more to say about this in §|| 

Tests 4 and 5 above drive home the point that the bias of an estimate of the un-normalized power spectrum from 
continuum-fitted data is highly variable. It depends on the redshift, resolution as well as S/N of the data. 

There have been in the literature discussions of an alternative method to normalize the quasar count: nor malizing by 
the maximum value of the continuum-fitted count, instead of by the mean count (e.g. McDonald et al. 1999| ). Note that 



this procedure is also sensitive to the S/N and resolution of the data. For instance, max(NQ/N§) = 1.4 in Fig. |6|, while 

max(Aq/7Vg) = 1.12 in Fig. [|, where Nq is the estimated continuum - they share exactly the same underlying cosmic 
signal but the former has a higher level of discrete fluctuations and poorer resolution - for reference, the true ma ximum 



trans missio n should be 0.99. This means one should take care to simulate the noise properties correctly (e.g. Ranch 



et al. |1997|) 



Lastly, we should emphasize that while trend-removal seems to be more desirable than traditional continuum-fi tting 



for the particular application here, continuum-fitting is still very useful for other purposes, which we will discuss in §4.2.2 
and §|^. But a fully automated procedure for continuum-fitting is clearly desirable. 

4.2.2. A Bonus of Trend-removal - Power Correction on Large Scales 



As is clear from some of the previous tests in §4.2.1, the power spectrum measured on large scales (i.e. scales comparable 



to the typical scales where the continuum has fluctuations) could contain spurious contributions from the continu um, th e 



is 



size of which depends somewhat on the continuum/mean-shape- model one assumes. The strategy adopted in the § |4.2.1 
a conservative one: assume a model for the continuum that is as simple (or smooth) as possible, perform trend-removal, 
and the resulting power spectrum would reflect the true transmission power spectrum at least on small scales, but not 
necessarily on large scales. 

Can we do better? The answer is yes, under certain as sump tions which we will make explicit shortly, and it illustrates 



an added benefit of trend- removal as introduced in the £ 4.2.1 . Readers not interested in the details can skip directly to 
the end of this section where two examples of how our technique of power-correction works are given (Fig. ^ and [To|) . 

Let us start by recalling the power spectrum estimator in eq. ([12]), but focusing now on the fact that the true N a is 
unknown, and has to be estimated using eq. (p"8|), which assumes implicitly that the true mean count obeys eq. (fRf), 

which is of course only a reasonable guess. Let us use Nq to denote the estimated mean count, which generally differs 

from the true mean count Nq. We have used N a somewhat sloppily before when we really meant Nq e.g. eq. |18[ ]. In 
other words, eq. |l8) should be more accurately written as 

Nq =J2 La0 K < L = P T (PP T ) _1 P (19) 



where p represents the basis functions. Similarly, the estimator for the power spectrum in eq. ( |L2] ) should be more 
accurately written as: 

£,(«) = J2 w afi {u)djdj , P 2 {k) = J2 w aP {k)d a f dP - b(k) (20a) 

a, P a,f3 

dj ee (A^ - N§)/N§ = D^N^/ ^ L^Nq 1 , D a '< = 5 a ~< - L a ~< (20b) 



We now assume the following quantities are small: dj and (^^ a7 -^Q — Nq)/Nq. The second quantity tells us how 
far off our estimate of the mean is from the true mean, while the first c ontai ns contrib ution s both from the fluctuation 
in transmission and from the second quantity. Therefore, putting eq. ( 20b ) and eq. ( [20a| ) together, the lowest order 
contributions to the expectation value of the estimator P2(k) are: 



(A(*)> 



V w^ikYD^D^ _ Q ^ )(1 



a/3777 

Pc(k) 



Q IV Q 



dk' 



P(k')G n (k,k r ) 



P c (k) ee V w af, (k)(D ar *D fh > _ Q J% ) 



(21a) 



(21b) 
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G„(k,k') = w aP {k)e lk ' ^ - uV) {D a ~< D 



otft-fq 



), 



(21c) 
(21d) 



where we have retained the old definition of 8J as (Nq 



(eq. 



2| ) . The above gives an idea of how biased the 

estimator P2(k) is. Note that we have used ( } here to include, in addition to the ensemble averaging as explained in §||, an 
averaging over the ensemble of possible continua (which changes Nq because it is directly proportional to the continuum 
count Ng). We have assumed the fluctuations in the continuum are uncorrelated with fluctuations in the cosmic signal 
SJ. We have also ignored the shot-noise contributions (e.g. b(k)) and will continue to do so for the rest of this section, 
because the scales where the continuum contamination could be a problem are typically large enough that shot-noise is 
subdominant. 

The term Pc{k) can be viewed as the power spectrum of the continuum fluctuation. This is fluctuation in the sense of 
D^Nq — Nq—L^Nq. This fluctuation would vanish if our trend-removal procedure were so accurate that the continuum 
shape is exactly reproduced. The term G n (k, k') is the effective window function, replacing the one in eq. (|l5|), which does 
not take into account the error involved in trend-removal. The desirable normalization condition J G n (k, Wjdk' / (2tt) = 1 
no longer holds with the choice of w a/3 in eq. (14a). We have instead 



/ G ^ fc ')£ = jf E R a *m^^ , 1 + e G (k) 

a/37 Q Q 



(22) 



where R a ^(k) is defined in eq. (14a). 

Assuming for now Pc{k) and tc(k) can be measured from a suitable ensemble of continua, we propose the following 
alternative estimator to P%(k), which removes the bias due to continuum contamination: 



P 3 (k) = [P 2 (k)-P c (k)]/[l + e G (k)} 



(23) 



The above gives us an unbiased estimate of the windowed power spectrum. The window is effectively G n (k, k')/[l+eG(k)] 
which has the desirable normalization. We will not attempt further improvements such as deconvolution in this paper. 
A useful alternative estimator, in cases where Pq dominates the bias in P2(k), is 



Pi(k) = P 2 (k) - P c (k) 



(24) 



The above estimator gets rid of most of the bias in the estimator P2(k) if Pc(k)/ P(k) ^ ec(k)- An interesting corollary is 
that, under such a condition, the bias in P2(k) is positive since Pc(k) is positive definite. Needless to say, this statement 
breaks down if Pc{k) is not the dominant source of bias, or if the fractional error in the mean-estimation is large (see e.g. 
Fig. |). 

It is intere sting to compare our derivation abo v e with the well - known one for the in tegral constraint bias in galaxy 
surveys (s.g. Peebles 198C; Landy & Szalay 1993; Bernstein 1994; Tegmark ct al. 1998). The integral constraint arises 
from the fact that the mean density of a galaxy survey has to be estimated from the same survey from which one is also 
measuring the power spectrum. The fact that the power spec trum estimator involves a non-trivial nonlinear combination 
of the data gives rise to a bias (see |Hui fc Gaztanaga 1999] ), quite analogous to our derivation here. However, in the 
standard derivations, it is assumed the shape of the mean density is known (often taken to be uniform), and therefore 
Pc(k) effectively vanishes, whereas ec(k) can be no n-neg ligible on scales compar able t o the size of the survey, but is 
otherwise small. The reader is referred to Bernstein (1994) and Hui & Gaztanaga (1999) for discussions on higher order 
contributions to the integral constraint. 

How should one estimate Pc(k) and ec(k)l Given an ensemble of continua (wit h cou nts represented by Nq), our 
procedure is to replace Nq, which appears in the definitions of Pc(k) and ec(k) (eq. plhfl , J22|) with iVg, and compute 
the corresponding ensemble averages. Note that Nq = N^(e^ T ), but (e~ r ), which is taken to be constant over the finite 
redshift range of interest, is divided out in the relevant definitions of Pc(k) and ec(k). 

The hard question is of course how to obtain a suitable ensemble of continua. The first thing one might try is to measure 
the power spectrum of the continuum-fits (i.e. Pc(k), or more generally, both Pc(k) and ec(fc)) from exactly the same 
regions from which one attempts to measure the transmission power spectrum. While this can give us a crude idea of 
how significant the continuum power spectrum is, it is not entirely satisfactory because part of what has been ascribed to 
the continuum might actually be large scale fluctuations in the cosmic signal that we are after, or vice versa. 

The second option that comes to mind is to measure the continuum power spectrum from regions where the continuum 
determination is relatively secure. Two possibilities are A. low redshift quasar spectra where the forest is much less 
dense, and B. regions of spectra which lie redward of Lya emission. The working hypothesis is that the continuum power 
spectrum in these two regions is the same as, or at least similar to, the one in the region where we attempt to estimate 
the transmission power spectrum (the forest of interest). There is no guarantee that the hypothesis is valid. For instance, 
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regarding possibility A, the continuum power could systematically evolve with redshift. In fact it probably does: assuming 
that the statistical properties of the quasar continuum in rest frame do not evolve with redshift, the observed continuum 
power would evolve as Pc[k, z\\ = Pc[k(l + zo)/ (1 + zi), zq]. One could in principle constrain such redshift evolution with 
a sufficiently large sample of low redshift quasar spectra. Regarding possibility B, it is not unreasonable to expect that 
the continuum power is higher on the red side c ompared to the blue side, because th ere are generally more broad emission 



lines on the red side (see e.g. Peterson 1997; Blandford, Nctzcr & Woltjer 1990; see below for caveats and a counter 



example, however). An upper bound on the blue continuum power is by itself interesting because one can then obtain a 
conservative estimate of how much spurious power is introduced by the continuum into one's forest power measurements. 
Furthermore, systematic differences between the red and blue continuum power can be studied and quantified with a 
sufficiently large sample of low redshift quasars. 

In Fig. ||, we show the continuum power spectrum measured from the continuum estimates on both sides of the Lya 
emission of a quasar at z=3 (QSO 1157+3143). The continuum estimates are shown in the bottom two panels. After 
fitting a flat mean to each echelle order, we compute the continuum power spectrum just as if this were the forest, and 
the results from the red side and blue side are shown as solid and dotted lines respectively in the top panel. The two 
power spectra look similar. However, we emphasize that because of the lack of small scale power in the continuum, most 
of the power on small scales (fc > lA _1 ) that we see in Fig. |s] is likely aliased from large scales. We will not attempt to 
perform a deconvolution to obtain the true small scale power; it suffices to note here that the true small scale power can 
only be smaller than what is shown in the figure. Also shown in the top panel is ec(fc), on both sides of Lya, which are 
basically indistinguishable from each other. Note that eo(k) <C 1. The second panel from the top shows the fractional 
difference between Pq from the red and blue sides, which is about 10%, with the blue continuum power systematically 
higher than the red one. The results here, though drawn from admittedly a very small sample, are quite interesting for 
several reasons. 



• The excess of the blue continuum power spectrum over the red one is consistent with the following hypothesis: that 
some of the fluctuations in the forest have been wrongly assigned to the continuum during the continuum-fitting 
process on the blue side. In other words, the true blue continuum power spectrum should be lower than the dotted 
line in the top panel of Fig. ||. An upper bound on the true blue continuum power spectrum is already very 
useful. One can use it to quantify how much, and on what scales, one should worry about spurious continuum power 
introduced into estimation of the transmission power. One can compare Fig. ^ with the theoretical expectation in 
Fig. [l], and see that the spurious power must be negligible for fc > 0.3A -1 . This explains why the determination of 
the transmission power spectrum from both the continuum-fitted data or the trend-removed data is very accurate 
in examples like Fig. ||, as long as one considers the normalized power. Unfortunately, the pieces of continuum 
we examine are not long enough to yield useful information on larger scales or smaller fc's. If one takes a crude 
extrapolation, the continuum power spectrum (or more accurately its upper bound) might become non-negligible 
compared to the transmission power spectrum at k ~ O.lA -1 . However, one must keep in mind that the theoretical 
transmission power spectrum in Fig. u\ is likely underestimated at small fc's because the simulation lacks large scale 
power. Nonetheless, there should be a genuine flattening of the transmission power spectrum at large scales. In 
any case, the first point to bear in mind is that an upper bound on the continuum power spectrum is useful as a 
conservative estimate of the possible spurious power, [j 

• Further, one can test the hypothesis that the excess in blue continuum power is due to contamination from the 
forest: if this is true, one expects the red and blue continuum power spectra to converge, as one goes to lower 
redshift quasars, because presumably, the blue continuum power spectrum should be less affected by the forest at 
lower redshifts. Even if their difference does not converge to zero (as suggested by the larger number of broad 
emission lines on the red side), but to some small but finite value, this is still a useful exercise because it gives us 
an idea of how different the red and blue continuum power spectra can be. If we can determine the blue continuum 
power spectrum to an accuracy of 10% say, and use this to correct for the transmission power spectrum on large 
scales, this is already a significant improvement over not correcting for the large scale power, or simply throwing 
away the information on large scales altogether. For instance, if the blue continuum power does become comparable 
with the transmission power at k ~ O.lA -1 , not subtracting off the spurious power would result in a fractional error 
of 100%, while subtracting off an approximate blue continuum power accurate to 10% reduces the error by an order 
of magnitude. 

Obviously, more testing using observed data is warranted, particularly on the estimation of red and blue continuum 
power as a function of redshift. This will be carried out in a separate paper. One natural question that might occur 
to the reader is whether a universal continuum power spectrum actually exists, given the large observed variations in 
the continuum from one quasar to another. It suffices to note that given an ensemble, the averaged power spectrum is 
always a well-defined quantity. The tricky part is to make sure the ensemble from which one estimates the continuum 

12 One should be aware of a possible pitfall of the above argument, however. It is not impossible that the opposite can happen, that one 
underfits the blue continuum, and ends up underestimating the blue continuum power. This is probably not the case here, where the data 
from which the blue continuum is estimated have high signal-to-noise and resolution (similar to the simulated spectrum in Fig. W, where it can 
be seen that the continuum fit tends to have features that follow the forest). Underfitting the blue continuum is more likely forlow resolution 
data, although even there, the situation is not clear: underfitting would result in underestimation of the continuum power on small scales, but 
not necessarily on large scales. Obviously, more tests are needed. 
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power spectrum has the same averaged continuum power spectrum as the ensemble of continua in the forest regions of 
interest. As a simple example: one might want to make sure the same proportion of radio-loud quasars are included in 
both ensembles. This is probably desirable if one uses the working hypothesis that low redshift blue continuum power 
is similar to high redshift blue continuum power, as suggested above. Alternatively, if the hypothesis that blue and red 
continuum power spectra resemble each other irrespective of redshift turns out to be a reasonable one, the simplest way 
to make sure one has the right ensemble is to use both sides of Lya for any give quasar: use the blue side for its forest, 
and the red side for its continuum. 

With all of the above caveats in mind, let us illustrate the technique of power correction with two idealized examples, 
where it is assumed the right ensemble of continua is in hand. 

In Fig. |^, we show in the bottom panel a simulated spectrum with a somewhat unusual continuum (middle panel) with a 
fair amount of fluctuations. We generate a set of 10 different continua and impose each on our underlying cosmic absorption 
to obtain a set of 10 different simulated spectra (only one of which is shown in the figure). We compute the power spectrum 
using P2(k) as in eq. (p0[). The resulting fractional error from the true transmission power spectrum is shown as a solid 
line in the top panel. There is clearly a lot of spurious power on large scales due to the imperfectly estimated mean count, 
which reflects the wild fluctuations in the continuum. We then apply the power spectrum corrections: the dotted line 
shows P3(k) from eq. ( p3|) while the dashed line shows Pi(k) from eq. (Q). One can see that subtracting the continuum 
power spectrum Pc(k) alone removes most of the spurious power. 

To make the example realistic, we have multipled the continua in the forest region by a power law that goes like (u a )°' 9e 
(i.e. the 'blue' continuum), and similarly multiplied the continua from which we actually estimate the continuum power 
by a power law of (u Q )~°' 01 (i.e. the hy pothetical 'red' co ntinuum). This is meant to mimic a possible turn-over of the 
quasar continuum around Lya (see e.g. ^hcng ct al. 1998] for evidence of a turn-over around Ly/3). We have in mind a 
situation in which the continuum power spectrum Pq is estimated from the red side of Lya. Clearly, the fact that the 
mean trends on the blue and red arc different docs not present an obstacle. 

In Fig. [l^ we show a similar version of the above, but with much noisier data and poorer resolution, and a mean 
power-law of (u Q ) 15 and a mean of (w Q )~ ' 9 have been imposed on the continua on the blue and red sides respectively. 
The same technique works here as well. 

One last point we should make: when the quantities Pc{k) and tc{k) are estimated from some ensemble of continua, 
they in general receive shot-noise contributions. We have ignored shot-noise here, assuming the scales where power 
correction is most interesting are sufficiently large that shot-noise is unimportant. This should be checked on a case by 
case basis. 



4.2.3. Gaps and Metal Absorption Lines 

There are at least two other possible sources of systematic errors in addition to that due to continuum-fitting. Gaps 
are quite common in observed spectra due to defects in the CCD, incomplete spectral coverage, or cosmic ray hits. 
Fortunately, since they are at known locations, we can either consider only those parts of the spectrum that are between 
the gaps (for instance when the gaps are large), or interpolate to fill in the g aps (f or instance when the gaps are small). 



The latter is what we have implicitly done in all of the tests mentioned in §4.2.1, where 3% of the pixels are assumed 
discarded because of cosmic-ray hits. The hits are typically one to a few pixels wide, and we simply fill them in by linearly 
interpolating the counts from neighboring pixels. Clearly, we can recover the power spectrum quite well in spite of the 
need to interpolate. 

A more challenging problem is possible systematics due to the presence of metal absorption lines. Shown in the bottom 
panel of Fig. |ll| is a simulated spectrum with resolution and S/N very similar to that of Fig. except that metal 
absorption lines as shown in the panel above have been added on top of the Lya forest. This list of lines is taken from 
a quasar spectrum redward of Lya (HS 1103+6416, z = 2.19). The mean transmission is estimated by assuming a flat 
trend for each echelle order as before. What is interesting is the dotted line in the top panel, demonstrating the creation 
of spurious power by the metal lines. The dashed / solid line shows fractional error in the power spectrum estimate if 
all metal lines with r > 0.4 / r > 1 are assumed "detected", and therefore cut-out and treated as gaps as before (i.e. 
interpolated across). Such a procedure seems to eliminate much of the spurious power. In practice, sufficiently strong 
metal lines should be identifiable by their narrow widths. 

Fig. O shows that metal absorption lines in data with lower resolution and poorer signal to noise have a relatively 
small effect on the accuracy of the power spectrum estimation. 



4.3. Random Errors: Shot-noise-bias, Variance and Minimum Variance Weighting 

Random errors arise firstly from (cosmic) sampling fluctuations, and secondly from electron /pho ton counting, which 
can be traced to fluctuations in the intrinsic quasar counts, the sky counts and the readout (see §[OJ). We will summarily 
refer to the latter as shot-noise. Shot-noise affects two aspects of power spectrum estimation. 



First, shot-noise introduces a bias which has to be subtracted off. This is the term b{k) in eq. (14b). We will give 
here a more general expression for b(k) suitable for different weightings {w a ^). As we have demonstrated in Fig. |j| 
shot-noise-bias subtraction can be important for low S/N data. We will return to this point in §|j| 

Second, shot-noise, together with cosmic fluctuations, determines the variance of the power spectrum estimate. We 
will give the expression for the variance in this section, and then address the question of how to best combine data with 
different levels of S/N to minimize the variance. 
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The power spectrum estimator we will focus on is given in eq. (|12[). It is assumed trend removal as explained in 
§4.2.1 has been performed. We ignore uncertainties due to the unknown continuum in this section. Here we do not limit 
ourselves to the choice of uniform weighting (eq. jL4j) as we have done so far. Let us rewrite w a ^(k) in eq. (|j~2| ) as 



lo^jfe) = w af} (k)R a0 (k) 



(25) 



where R a ^{k) is given in eq. (14) and is an average of the Fourier basis over some bin centered at k, with width Afc. 
It can be shown that the variance of such a bin-averaged power estimate is given by (Appendix B) 



v(k) = ([P 2 (k) - p(fc)] 2 ) = £K«(fc)] 2 £«(fc) 
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(26b) 



N a 
Q 

n a 



c_ 



assuming the k of interest satisfies k 3> l/C and that the width of the bin Afc also satisfies Afc ^> 1/C, where C is 
the length of the spectrum. This i s sometimes refe rred to as the classical limit in the case of galaxy power spectrum 
measurement (Feldman ct al. 1994; Hamilton 1997a). We will not consider larger scales here, because measurements on 
such scales are likely dominated by systematic rather than random errors. 

The symbol n k denotes the number of m odes within the k bin of interest, N is the number of pixels in the length C; 
q a , iVq and V£ are as defined in eq. (|l4b|) . The quantities q' a and q" a are analogous to q a : 



q 



ict a 

ps 9b 



r = E 



(27) 



The symbols (T)k k , (B)kk an d {P)kk represent the shell-averaged trispectrum, bispectrum and power spectrum respectively 
(( ) here is to be distinguished from ensemble average discussed in §2|): 
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(28a) 



(28b) 



(28c) 



where the sum over k[ extends over modes within the bin centered at fci, and similarly for fc 2 . The trispectrum T and 
bispectrum B are Fourier transforms of the four and three point correlation functions, defined in an analogous manner to 

eq " @- ■ 

The variance as given in eq. (Eq) contains contributions from both cosmic fluctuations and discrete fluctuations (see 

§||). The terms such as P(fc) 2 and (T)kk arise because of intrinsic fluctuations of the cosmic signal from one part of 

the universe to another - these terms are present even if one has data with arbitrarily high signal-to-noise. The terms 

containing Nq arise because of discrete fluctuations - these we will loosely referred to as shot-noise. 

As we have emphasized in §tO & 4.1, the shot noise contributions to the random erro r are not exactly Poisson- 

distributed. The shot noise contributions (ignoring cosmic sample fluctuations) in eq. ( 26b ) would all be simply 1/Nq 

if TVq were strictly a Poisson variable. We have additional fluctuations in Nq due to the background (sky and readout), 
and also due to non-unity weights used in reducing the data (eq. ||; see also end of §3.1). 

Given eq. (|2"6|), it is simple to derive a weighting w a ^{k) that minimizes the variance V{k), subject to the constraint 
that the effective window (G as defined in eq. p5|) is properly normalized. This is most simply derived by minimizing 
the following Lagrangian: 

L(k) = V(k)-X(l G(fc,fc')^-1) (29) 
where A is a Lagrange multipler. Differentiating the above respect to w a ^(k), and setting the result to zero, we obtain: 

w af3 {k) = [P Q (fc)£; /3 (fc)]-V M ( fc ) , M(k)= y ^[E> 1 (k)}- 1 M/C (30) 
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The corresponding shot-noise subtraction, instead of eq. (|14b|), would then be 



b(k) 



\k) 



q a N§ 



(31) 



where q a , Nq and Vg are as defined in eq. (14b). 

In summary, the minimum variance estimator of the power spectrum is 

h(k) = ^2R aP {k)[E a {k)^5J][E {k)^5 l3 s ]/M{k) - b(k) 

a, (3 



(32) 



where b(k) is given by eq. (|3l"|), E a {k) and M(k) are given in eq. ( p6|) and (|30|), and R al3 {k) is as in eq. (14a). 

The minimum variance estimator can be understood simply as follows: 5^ at each pixel is weig hed by l/y/E a {k) before 
the array is Fourier transformed, squared and grouped to form band power estimates. Note that the above estimator 
reduces to the one with uniform weig hting (eq. |1) if E a (k) were independent of a, e.g. when sample/cosmic variance 
is significantly larger than shot-noise (P(k) [CjJ\][q a Nq + Vq]/[Nq] 2 ). It is important to note that the weighting as a 

function of a is determined by Nq rather than, say Nq. Down- weighing pixels with a lot of absorption (hence relatively 
low Nq) would be a wrong thing to do, since the fluctuations in absorption is the signal that we are after. The proper 
procedure is to down- weigh pixels with an overall lower mean count Nq. 

Unfortunately, the minimum variance weighting given above is difficult to implement because one needs to specify 
simultaneously P, B and T, in addition to the level of shot-noise. A common simplification is to use the Gaussian 
approximation in which E a {k) is approximated as: 



E a (k) 



P(k) 



C q a N§ + V£ 



n (N°y 



(33) 



(see e.g. Hamilton 1997a). Note that in addition to ignoring the bispectrum and trispectrum terms, the above also ignores 
certain power spectrum terms which are mixed with shot-noise - the last three terms in eq. (|26|), which is equivalent 
to assuming that either the shot-noise or the correlation is sufficiently weak. With the above approximation, one can 
start with so me initial P and us e the minimum variance weighting scheme to get a first measurement of P, and iterate 
subsequently (Bond et al. 1997). Analogous (Gaussian) power spectrum estimators for galaxy-surveys and microwave 



background experiments have been w i dely discussed in the literature (e.g. Feldman et al. 1994 ; Vogeley & Szalay 1996 
Tegmark ct al. 1997| |Hamilton 19"97a]; |Tegmark et al. 1998|; |Bond et al. 1997t |Scljak 199. 



We will not attempt to address here the important question of how significant the non-Gaussian contributions are. A 
proper treatment will involve the analysis of a large number of simulations or a large data-set, which we hope to present 
in a future paper. It suffices to say that the very nonlinear mapping from the density field to e _T will likely introduce 
some degree of non-Gaussianity, even if the initial density field is Gaussian. 

The use of observed data to study this issue is particularly interesting, and deserves some more comments. In principle, 
since different QSO sightlines typically sample independent regions of the universe, one can estimate the variance of 
the transmission power spectrum, and hence infer the importance of the non-Gaussian contributions, using directly the 
fluctuations in power spectrum estimates from one sightline to another. However, one should keep in mind that shot- 
noise also contributes to the variance. Since different lines of sight generally have different S/N, the sightline-to-sightlinc 
fluctuations in power spectrum estimates should be interpreted with care. In a dataset of several quasars, it is possible that 
the quasar-to-quasar fluctuations are dominated by a few with low S/N, and their mean-square would give an overestimate 
of the true power spectrum variance. 

We show in Fig. [l3| an example in which the data consist of 6 high quality spectra (similar to Fig. ||) and 6 others 
with S/N about 20 times smaller. The bottom panel shows the power spectrum estimated with uniform weighting (eq. 
fl4|| ) while the top panel represents the power spectrum estimated with minimum va rianc e wei ghti ng using the Gaussian 
approximation. The (1 a) error-bars are theoretical - they are estimated using eq. (26a) and (|33|). This illustrates how 
our weighting scheme can reduce the error bars at high fc's where shot- noise is important. Absent information on the 
non-Gaussian nature of the power spectrum variance, we advocate the Gaussian weighting scheme (eq. p3| ) as a rational 
way to combine data with different levels of S/N to reduce the variance, even though it does not necessarily achieve 
minimum variance. 

In combining the different spectra with different S/N, we have weighed the power spectrum estimate of each line of 
sight by its inverse variance, which is an obvious generalization of the minimum variance weighting introduced above. For 
instance, suppose we have two separate lines of sight A and P, we could combine the two power spectrum estimate P£ 
and PP in the following way, assuming the two lines of sight are independent: 



A(fc) 



P£{k) 



V A (k) V B (k) 



V A (k) V B (k) 



(34) 
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where V^fc) and Ve(fc) are estimated with the same P but could have different levels of shot- noise. The noisier quasar 
spectrum is naturally down-weighed. 

Lastly, we should emphasize the above discussion does not address the issue of cro ss- variance between power spectrum 



estim ates a t two differen t wave bands, which is introduced by the non-Gaussian terms (Mciksin & White 1999; Scoccimarro 
1999] ). Hamilton (199E) introduced a scheme which simultaneously diagonalizes the covariancc and minimizes it. 



et al. 

However, it makes specific assumptions about the form of the trispectrum and bispectrum, the validity of which for the 
forest remains to be checked. 



5. DISCUSSION 

Our recipe for measuring the transmission power spectrum is summarized here. 

• Given an array of reduced quasar counts Nq, identified metal lines should be removed, especially the strong ones 
(t > 1). Small gaps in the spec trum (e.g. due to cosmic-ray-hits removal) can be (linearly) interpolated over, while 
large gaps should be avoided (§ 4.2.3 ). 



The mean quasar counts (Nq) is estimated using eq. (|18[). The mean-basis (the functional form of the mean-trend) 
should be chosen to be as smooth as possible - underfitting is better than overfitting (see Test 1 and 2 of §4.2.1). 
In practice it appears that a flat mean suffices for short spectra (~ 50 A), while polynomials up to the third order 
can be used for longer spectra (~ 500 A). One can gain an idea of what a reasonable basis is using the red side of 
Lya or low redshift QSO spectra where the continuum can be seen more clearly. 



• Define the trend-removed and normalized fluctuation 5" according to eq. (12b), and the power spectrum is estimated 
using the quadratic estimator in eq. ( ]12a| ). Different weightings are possible, the simplest that we recommend is 
given in eq. (1 4| ). A more sophisticated weighting scheme which can reduce the random error is given by eq. (p0|), 
( |3l| ) , ( j3 2| ) and (p3|) . If one is interested in the real-space correlation function instead, the recommended weighting is 
eq. (|13|) - this gives a smaller varianc e com pared to other estimators commonly used in the literature. We emphasize 
that the shot- noise-bias (b(k) in eq. ]14b[ , or in eq. pl]|) should be subtracted correctly, especially for noisy data. 

• If better control over syste matic errors on large scales ( ^ 30A) introduced by the unknown continuum is desired, 
the techniques outlined in § 4.2.2] can be used. The c orres ponding estimator is given in eq. (|23|), which requires an 
estimate of the continuum power spectrum Pq (eq. |21b|| ) and an additional correction fact or eg (eq. p2| ). This 
procedure requires the identification of an appropriate set of continua (see discussions in S 4.2.2 ) . Even if one is 
not interested in the power spectrum on large scales, we recommend this procedure as a consistency check that the 
spurious power introduced by the continuum is negligible on the scales of interest. 

What implications does the above have for one's observing strategy? To discuss this question, we need to take a closer 
look at the issue of shot-noise. The shot-noise enters in two different places in the above discussion. First, it contributes 
to the variance (random error) of the power spectrum estimate (eq . [|26|). Second, it appears as a bias in the power 
spectrum estimate that we have to subtract off (e.g. b(k) in eq. [12a | fc |14[|). 

In the lit erat ure on power spectru m measurement, shot-noise subtraction has been largely ignored (e.g. |Croft 



et al. 1998a ; see McDonald et al. 199E for an alternative approach where shot-noise is simulated rather than subtracted). 
Let us estimate how important it is. The expression in the simplest case of uniform weighting is given in eq. ( |l4b| ) (see 
eq. f3l|| for more complicated weightings), which can be rewritten as 



.... _ Au ^ E,(^ a ) 2 [ff^ffg^g 



(35) 



where Au is the size of a pixel, TV is the total number of pixels, and the rest of the symbols are as defined in §3.1 



i is 



the pixel-label in the spatial direction and a in the spectral direction, Nq is the mean reduced quasar count, V^ a is the 
background variance, W la is a weighting, and g^ and g£ are the point-spread function and blaze function respectively 
(eq. H & Q). An important observation is that the numerator within the summation is closely related to the variance 
array which is often given along with a spectrum (eq. 



v a r(a) = ^(Wn 2 [N^+V^} 



(36) 



The quantity Nq is of course different from g^g^Ng which we need to estimate the shot-noise, but since we are in 
practice interested in an average over all pixels, it turns out the following estimate of the shot-noise is accurate to within 
a percent for all cases we have tested: 

var(a) 



3 See paragraph after eq. (H) and [Hon the distinction between N^* and N t 



Q 



16 



Without the above approximation, an exact estimate of the shot-noise would require the knowledge of Nq, J2i(^ a ) 9ps9b 
and ^(W^) 2 ^". 

Eq. (|37| ) provides a useful means of estimating the shot-noise (see Appendix A on shot-noise estimation under more 
complicated circumstances i.e. with non-trivial rcbinning or weighting). One can simplify further by making a crude 
approximation in relating b(k) to the typical signal-to-noise ratio of the data (which is often quoted at the continuum) 
through 

6(A) ~ (Au/f)(N/S) 2 (38) 

which can be justified if one ignores the part of the variance due to the sky and readout. We find that this simple rule 
of thumb generally provides an underestimate of the shot-noise (particularly at low S/N where the background counts 
become important), but is accurate to within about a factor of 2. 

Fig. summarizes some useful information for devising an observing strategy, based on our estimate of the shot-noise 
in eq. (|3q) a bove. The solid line shows the mean observed transmission power spectrum at z = 3 taken from McDonald 
et al. ( |1999| ).p| The two horizontal dotted lines show the level of shot-noise expected for the 2 extremes of the kinds of 
observations we are likely to e ncounter - the bottom corresponds to v e ry high signal to noise (S/N) observations with 



HIRES quality resolution (e.g. Hu et al. 1995; Kirkman fc Tytler 1997; Rauch et al. 1997), while the dotted line on the 



top corresponds to low S/N observations expected for a large number of quasars in the Sloan Digital Sky Survey (SDSS). 
We emphasize that the shot-noise level does not depend on the resolution per se, but on the pixel size for a given S/N. 
SDSS is expected to produce ~ 1000 QSO spectra at S/N = 20 per pixel (QSO's at z > 2.7, where the redshift limit is 
determined by the blue limit of the spectrograph, 38 00A), ~ 10000 at S/N = 15 and ~ 30000 at S/N = 7, corresponding 
to i' magnitude-cuts at 18, 19 and 20 respectively ( Fukugita et al. 1996| ; Fan 1999] ). The pixel size of SDSS is quite 



uniform in velocity 70km/s, which is equivalent to 1.133 at 4864A (Lya at z = 3). 

Clearly, the importance of shot-noise depends on the scales at which one is interested in measuring the power spectrum. 
A few interesting scales are shown in Fig. [l4|. First, instrumental resolution imposes a high k limit beyond which one 
cannot reliably measure the transmission power spectrum. The resolution window is often characterized by a Gaussian with 
a given FWHM. The effect of such a resolution window on the power spectrum can be represented by P(k) — > P(k)e~ k l k ~° 
where k a = y/2/h\2/ FWHM ~ 1.7/ FWHM . Two representative k a 's are shown as long tickmarks at the top. Note that 
even at k = k a , the resolution window reduces the power by 63% and so has a non-negligible effect. The Sloan FWHM is 
about 2.1 pixels i.e. 147km/s, or 2.4A at 48641. 

On the other hand, the range of scales that is currently being used to infer the mass power spectrum is indicated by 
the interval near the bottom Afc usc . The high k limit is set by the scales at whi ch the shape of th e power spectrum is 



preserved in the transformation from mass to transmission (i.e. linear biasing e.g. Croft et al. 1998a ). We can see that for 
high quality Keck spectra, information from a whole decade of measurable scales is unused for the recovery of the mass 
power spectrum - it would be very useful to push the current analysi s techniques to t hese scales, since power on these 



scales is of particular interest in constraining e.g. neutrino properties ( Hui et al. 1997 ; Croft et al. 1999| ). Such an effort 



would require disentangling the effects of peculiar velocities and thermal broadening, however. At the other end, the low 
k limit of currently usable scales is set by the scales at which the continuum fluctuates. This is indicated by the dashed 
line at the top, where the transmission power spectrum is unknown. 

From the above discussion, we can distill a few tips for observing/ analysis. 



To ensure that shot-noise is subdominant, one might want to achieve S/N > y 10Au///P(fci nt . ) where fci n t. is the 
scale of interest, and / is the mean transmission. The factor of 10 is somewhat arbitrary - this will ensure the shot- 
noise contribution to the power spectrum variance is no more than about 20% (under the Gaussian approximation; 
see eq. [|U), or tne 1°" error-bar on the power spectrum would only be increased by 10% due to the contribution 
from shot-noise. An important question is what fci nt . should be - that depends on at what scales one can usefully 
extract cosmologically interesting information. Current literature mainly focused on k ult . ^ 2A _1 , where P ~ 0.06A, 
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therefore S/N ~ 20 y Au/IA would be sufficient. Since P rises with scale, shot-noise would be even less important 

at smaller k's. Note that with very small Am such as ~ 0.05A, S/N as low as 4 — 5 is acceptable. To give some 
examples, a S/N of 8 per 0.05A can be achieved with an hour of exposure using Keck/HIRES for a V — 19 quasar; 
on the other hand, a S/N of about 15 per 1.1 A is expected with just slightly under an hour of exposure using the 
SDSS spectrograph for a i' = 19 quasar. 

• A corollary of focusing on only k < 2A~ 1 is that observations with k a > 3 x 2 a 1 or a resolution FWHM of 0.3A 
or R ~ 16000 at z — 3 are acceptable. The factor of 3 above (i.e. in 3 x 2A~ 1 ) is somewhat arbitrary - it ensures 
that at k = 2A _1 , the resolution window does not reduce the power by more than 10%. If the resolution window 
is known accurately, or if one is willing to sacrifice information on the small scales close to k ~ 2A , one could in 
principle consider lower resolutions. We would like to emphasize, however, that in principle, the modes at k > 2A~ X 
could still contain very interesting cosmological information, even though the current attempts at recovering the 
mass power spectrum ignored them. 

We divide the un-normalized power spectarrn (eq. jjj]) of Mr-Dona H et, al. by the square of their measured mean transmission to obtain 



the normalized power spectrum given in Fig. |14| See tests 4 and 5 in §1.2.] on the bias of the un-normalized power spectrum 
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If shot-noise is subdominant compared to the power spectrum, the only other limiting factor to the size of the 
random error is the total size of one's sample or the number of sightlines in it. Assuming all sightlines have 
similar coverage with length C, then the fractional error of a single k-mode (i.e. in a k-bin of 2it/£) is given by 
SP/P = C/^/Wsight where N si ^ is the number of sightlines assuming they are independent, and C — 1 under the 



Gaussian approximation (eq. [p3|), and a little larger than unity under more general circumstance (see e.g. Meiksin 
[| White 19991 ; fccoccimarro et al. 199S| ). 



• How should one distribute one's observing time among quasar targets to minimize the random error on the trans- 
mission power spectrum? There are many possible versions of this problem. We will discuss two, giving an explicit 
solution for the first, and only general expressions for the second. In the simplest case in which all the candidate 
quasar targets have similar magnitudes, given a finite amount of observing time, one can deduce the optimal total 
number of quasars one should target by 

minimizing iVtot." 1 [^(fcint.) + — ] 2 subject to iV to t.t = T tot . (39) 

/ t 

where iV to t. is the total number of quasars targeted, ki nt , is the scale of interest, X"tot. is the total amount of observing 
time one has, t is the amount of time one spends on each quasar, and 1/A is equal to (S/N) 2 reached per unit 
exposure time. The above assumes eq. ([33]) and that the sightlines are independent. The solution is easy to deduce: 
JVtot. = [P(h nt .)f/Au][T tot ./A], or t = AAu/f/P(k int .). A typical value for 1/A is 1/A ~ 1200hourper A x 
Iq(19~ mag.)/ 2. 5 x [ a p er ture/100 m. 2 ] x /throughput, where /throughput is about unity for Keck/HIRES, and ~ 2.5 for 
the SDSS. Using again k lnt . ~ 2A _1 , for a 19th magnitude quasar, with an aperture of 6.25 m. 2 and assuming 
Au = lA and /throughput = 3, the exposure time is t = 14 minutes only! The above prescription, however, only 
allows for just enough_ exposure time to reduce the shot-noise to a level comparable to the cosmic/sampling variance 
(i.e. P(fcmt.) ~ (Au/ f)(N/ S) 2 ) - the sole aim is to maximize the number of quasars observed within a given length 
of time to beat down the sampling variance. The prescription would certainly be different if one has, for instance, a 
finite number of quasar targets, or if one has other purposes in mind - such as measuring the mean decrement, etc 
(see earlier prescription for making shot-noise subdominant, equivalent to multiplying t by about a factor of 10). A 
more general version of the above problem deals with a case where the quasars span a range of magnitudes i.e. A 
is no longer the same number for each quasar. A simple ansatz is to assume t = aAAu/ f / P(k- nrt .) , in other words, 
spending more time for fainter quasars because it takes longer to beat down the shot-noise, except that we have a 
normalizing factor a which enforces the constraint of total observing time: a = T tot ,[fP(h m t,)/ Au]/ A'n(A')dA' 
where n(A)dA is the number of quasars with A falling in the given range, and A m i n corresponds to the brightest 
quasar in one's sample. Then, we can determine how many quasars one should include, starting from the brightest 
one, or how faint one should go by minimizing [J^ n(A')dA']~ l P(k- mt .) 2 (1 + 1/ce) 2 with respect to A. 

The following is particularly relevant for SDSS or comparable observations. 



In addition to contributing to the power spectrum variance, shot-noise also contributes a bias which has to be 
subtracted off (see e.g. Fig. ^). This is quite important for SDSS because, with > 10 4 sightlines, the survey has the 
capability of reducing the fractional error of the power spectrum to < 1% per mode. Therefore, a bias of ~ 3 — 100%, 
depending on the scale of interest (as indicated by the top dotte d line in Fig. |l4l), i s not acceptable and s hould be 



subtracted off. We note that analyses so far in the literature (e.g. Croft et al. 1998bj; McDonald et al. 1999) focused 
on higher quality data where S/N ~ 30, with Am ranging from about 0.04A to 1A, and so according to eq. (38) 
and Fig. n% the shot-noise bias was about 1% of the power or smaller and therefore could be ignored, although a 
more careful check should be performed for some datasets with lower S/N. 

The low resolution of SDSS spectra implies that it would be difficult to obtain useful information on scales fc 3> k a ~ 
0.01 (km/s) -1 or 0.7A -1 . On larger scales or smaller k's, two problems have to be reckoned with. For k ~ 0.3k a — k a , 
the resolution window suppresses the power by 10% or more - therefore, one needs to have an accurate measure of 
the resolution window to recover the true transmission power spectrum. This can be achieved by using narrow 
metal lines or arc lines. There are relatively few sky lines in the relevant part of the spectrum. 

For scales k < 0.0 04 (km/s) -1 or 0.2A -1 , the effect of the continuum has to be properly taken into account, and the 



method of §1.2.2 can be useful here. From Fig. u4, it is clear that the range of scales accessible to SDSS would be 



quite limited unless a correction for continuum contamination is applied. 

It is worth pointing out that while the above quoted numbers are all based on z = 3, we do not expect them to change 
gnificantly for z = 2 or z = 4. This is in part because of the slow evolution of the transmission power s pectrum 



le g rowth of the mass spectrum with time is partially compensated by the lowering of the mean decrement (McDonald 



al. 1999) 



3 We thank Rupert Croft for useful discussions on this point. 
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At least three issues remain to be explored in future work. As we have em phasiz ed, the concept of correcting for 



continuum contamination in the transmission power on large scales as laid out in §4.2.2 is an interesting one, but requires 
more testing. An important check is the measurement of continuum power spectrum as a function of redshift on both 
sides of Lya emission. Second, counts-in-cells analysis (i.e. measuring moments of the one-point probability distribution 
function), just like power spectrum analysis, requires shot-noise subtraction, and since typically one considers cells at 
the limit of resolution, shot-noise is likely non-negligible except for hig h S/N spectra. Counts-in- c clls analys i s provides 



a ver y interest ing way to test the gravitational instability paradigm ( |Gaztanaga fc Croft 1998|; Hui 1998 ; Nusser & 



Haehi icTt 1999] ), and should be done with care. Useful expressions will be presented elsewhere ( [Hui 2000| ). Lastly, as is 



hopefully clear by now, the spirit of the methods presented in this paper is to avoid continuum-fitting and replace it with 
trend-removal. We have demonstrated that this is possible for measuring the transmission power spectrum. However, 
other quantities of interest such as the mean decrement requires an estimate of the continuum, by definition. Furthermore, 
to theoretically interpret the transmission power spectrum in terms of the mass fluctuation, current methods require the 
measurement of the mean decrement to fix a free parameter in one's cosmological model, which is a combination of 
the ionizing background, the mean baryon density and the mean temperature. Therefore, in a sense, the technique of 
trend-removal only goes half-way in solving the problem of continuum-fitting. Although we still recommend our method 
over continuum-fitting because the transmission power spectrum is an unambiguous quantity which can and should be 
determined as accurately as possible (not to mention the fact that continuum- fitting is difficult with low S/N or low 
resolution, or at high z), there is clearly a need for an alternative method to bridge the gap between the me asured 



trans mission power and the theoret ically interesting mass power. This will be explored in future publications (Hui & 



Buries 2000; Zaldarriaga et al. 2000) 



We thank for useful discussions Lcn Cowie, David Kirkman, Patrick Petitjean, Michael Rauch, Wal Sargent, Don 
Schneider, David Weinberg, Don York and the participants of the 1999 Haifa workshop on the intergalactic medium and 
large scale structure. We also thank Nick Gnedin for supplying a simulation. Special thanks are due to Matias Zaldarriaga 
for pointing out the importance of aliasing and for many interesting discussions. This work was supported in part by the 
DOE and the NASA grant NAG 5-7092 at Fermilab, and by the NSF grant PHY-9513835. LH thanks the IAS for the 
Taplin Fellowship. 

APPENDIX A 

Our main aim here is to derive eq. ( |l5| ) for the estimator P2(k) which is given by eq. (|l^) and eq. (fl4|), with an eye 



towards generalization to W 1 ^ different from eq. (||) and W al3 ^ 8 a P . Derivations of results in §13 on estimator variance 
and the minimum variance power spectrum estimator can be found in Appendix B. We will ignore the integral constraint 
and assume N a is known to high accuracy. 

We need first of all the correlation matrix (8J5^). We will do it in 2 steps, first let us derive (5J5^)d- Rewriting 8f 



(cq. [12b|) as (N a - N a )/N a + (N a - N a )/N a where N a = {N a ) D , it can be shown that 



(Sffyo = SfS> + --p((^ - N a )(N? - N")) D (40) 
iv Q JV Q 

- SfS f } + ]T W^W^(W^) 2 ({N^) D + vp) 

Q Q 7* 

where 8J is to be distinguished from SJ in that it has only cosmic or sample fluctuations (eq. ||), Nq is a strictly Poisson 

variable with an average given by eq. (|5|), and V^ 7 is the variance contributions from the sky and readout (eq. |^). 
Taking the cosmic mean of the above, we obtain the correlation matrix 

= e(«° - «") + — ^ E w^w^ J2^ W ^ 2 &^ Q ) + ((w^vp] (4i) 

The second term on the right is the shot-noise contribution which has to be subtracted off. Using the estimator in eq. 
dl2]), with w af3 (k) given in eq. (|l4|), the correct shot-noise subtraction is: 

b(k) = {C/N 2 ) £ W^W^W^fg^slNl) + <(VF l7 ) 2 )^ 7 ]/(A^) (42) 
a/37 * 

where we have made use of the assumption that W ai is non-zero only for a and 7 on very small separations, and that 
the k of interest satisfies k(u a — u 1 ) < 1 on such separations. 

Note that for weightings such as the one given in eq. (Q), W n depends on Nq which makes an estimation of the shot- 
noise non-trivial. However, simplification results in two extreme limits: in the signal dominated regime where Nq 3> , 
W %1 reduces to uniform weighting as in eq. (||); in the background dominated regime where Nq Ug 7 , W %1 reduces to 

(l/3b)(Sps/^B 7 )/(Sj(5ps) 2 /^B 7 ) wn i cn is independent of Nq. In such cases, and for W ai = 5 ai 1 the above b(k) reduces 
to the one given in eq. (|l4b|) . 
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For W a ~< / 5 ai (i.e. rebinning has bee done; we continue to assume W %1 is roughly independent of the signal), we can 
write b(k) as 

b(k) ~ (C/^) £ ^^-^{W^flg^glNl) + V^] (43) 

a/3 7 iV Q iV Q i 

where the term under summation of i is simply the variance array of the pre-rebinned data (eq. p6| ) , and one can replace 
TVqTVq by (^Vq) 2 , since W^W^ 1 is non-zero only if a and j3 are close together. Eq. ( ^7| ) is therefore replaced by the 
following if the data have been rebinned: 

where var (a, /?) is the variance matrix of the rebinned data. 

To complete our derivation, we need to show that the choice of w a P (k) given in eq. (|14|) has the correct normalization 
such that the window function satisfies J dk'G(k,k')/2ir = 1 (eq. [jUSf). Putting eq. (|14|) into eq. (|l2]), and using the 
correlation matrix given in eq. (^T|) together with the relation between the two point function and the power spectrum 
in eq. (g), it is not hard to see that (P 2 (k)) satisfies eq. (fL|) with G(k,k') given by J2 a p w Q ^(fc)e lfc ' (tl =~ u ' 3) . Using 

J dk' e^ k _ 2n5 a Pjy / £ then comple tes the derivation. One might want to explore more complicated data windowing 

(e.g. Press et al. 1992b ; Hamilton 1997b ), but since in practice uncertainties in the large scale power estimate, where 
the survey window matters most, are likely dominated by the continuum, the simple choice we have adopted is probably 
adequate. 

APPENDIX B 



We derive here the band-power variance given in eq. (26). The power spectrum estimator is given in eq. ( [12] ) with the 
matrix w a ^{k) given by eq. (p5|). We ignore here the uncertainty in the estimation of the mean count Nq. 
The band-power covariance can be written compactly as 

C(k u k 2 ) ee {AP a (k 1 )AP^(k 2 )) = ]T w^(h)w^(k 2 )((5J5fS}5}) - {5f5 f }){5}5})) (45) 

a/3777 

The band-power variance is simply the diagonal piece of this matrix: C(fc, k). 

We can work out (S^SjSjS^j) — (SJS^^JS^) using the same methodology as used in Appendix A for (8J8^)\ rewrite 8J 

as {N a - N a )/N a + (N a ~ N a )/N a where N a = {N a ) D , and as before, take the discrete-ensemble average (} D before 
taking the cosmic average ( ) . The result is 

C(h,k 2 ) = ^ J (ki)w at3 (k 2 ) ( (46) 
a/371) \ 



£ w^w^((w i n 2 K a s 9^(i + 6?) + v^sjsf) 



+1 other perm. : (a <-> 7, (3 <-> 77) 

+3 other perm. :(/3 «-> a, 7 <-» r/), (7 <-> 77), ((3 «-> a) 

- £ ^^^r ff ((r ff ) 3 [ s »^^(i + ^) + A>r]^) 

iV Q A/ Q iV Q cri 

+3 other perm. :(a <-> /?), (a <— 7), (a rj) 

+ 1 £ ^^^^((^) 2 (iu«) 2 ( 3 - 5 ^5(i + ^) + vr) 

(#^(1 + 5*)+^*) 

+ 1 £ ^^^w^((iu-) 2 (iu«) 2 ( 3 - 5 -a7-(i + + v-) 

(diidb^Qi 1 + + ^b X ) + 1 otllcr perm. :(a <-> /3) 

+ , tt ,Ly, £ w^^tu" ct w^((^) 4 ^^^(i + «y)> 

JV Q iV Q iV Q JV Q ai 
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The first set of terms (second line) arise from the shot- noise- free part of (5j , namely (N a — N a )/N a . The next two sets 
of terms (third + fourth and fifth + sixth lines) come from combinations of 8J involving products of two shot-noise terms 
with two shot-noise-free terms. The next set of terms (seventh + eighth lines) arises from products of three shot-noise 
terms and a noise-free one. The next set of terms (nineth to thirteenth lines) comes from products of four shot-noise 
terms. The last term corresponds to what has to be subtracted off to compute the covariance. 

To make further progress, we assume W la is independent of 8J, which is strictly correct for W la given by eq. ([|), 
but only roughly so for eq. ((?]) (see discussion in Appendix A). Then, taking the small wavelength limit in the sense 
that k,Ak > 2tt/£ (Ak is the size of a £:-bin), and making use of the fact that W aa W^ a is only non-zero at separations 
u a — much less than 1/k where k is the wavenumber of interest, we obtain 



C(fc 1; fc 2 )=^(^) 2 ^^ 

^-i P (ki) + n^p((^) 2 (5p^b^Q + Vb)) E W^W^] 2 6 k ^ + 

+^(B) klk2 E 7^<(^) 2 «} E w^w^ 

+2(P(k 1 ) + P(k 2 ))^ E 7^p(( W/m ) 3 «) E W^W^W^ 



(47) 



011 



2 JP 



A few comments are in order. The terms in the third + fourth 

fourth lines of eq. 



nineth + tenth lines of eq. (|46|) are canceled by the 
) contain contributions proportional to B{k, — fc,0) 



last term of eq. (46). The terms in the third 

which vanish. There is also a term from the nineth line of eq. ([46) that is proportional to P (0) which vanishes also. The 
shot-noise terms in the second line of eq. (E7j) come respectivelyfrom the fifth + sixth + nineth + tenth lines of eq. (|46| 
The rest of the terms in eq. d47j ) basically follow the order t hey are presented in eq. (f46|). 
Lastly, setting W a/3 = 6 af3 and k\ = ki then recovers eq. (|2(j). 

REFERENCES 

Barlow, T. A. & Sargent, W. L. W., 1997, AJ 113, 136+ 

Blandford, R. D., Netzer, H., & Woltjer, L., 1990, Active galactic nuclei, Springer- Vcrlag 
Bernstein, G. M., 1994, ApJ 424, 569 
Bi, H. & Davidsen, A. F., 1997, ApJ 479, 523+ 

Bi, H. G., Boerner, G., & Chu, Y., 1992, A&A 266, 1 

Bond, J. R., Jaffe, A. H., & Knox, L., 1997, preprint |astro-ph 9708203| 



Bond, J. R. & Wadsley, J. W., 1997, Proc. 12th Kingston Conf., Halifax, |astro-ph 9703125 
Cen, R., Miralda-Escude, J., Ostriker, J. P., & Rauch, M., 1994, ApJ 437, L9 
Cen, R., Phelps, S., Miralda-Escude, J., & Ostriker, J. P., 1998, ApJ 496, 577+ 

Cristiani, S., D'Odorico, S., D'Odorico, V., Fontana, A., Giallongo, E., & Savaglio, S., 1997, MNRAS 285, 209 
Croft, R. A. C, Hu, W., & Dave, R., 1999, Phy. Rev. Lett. 83, 1092 
Croft, R. A. C, Weinberg, D. H., Katz, N., & Hernquist, L., 1997, ApJ 488, 532+ 

Katz, N., & Hernquist, L., 1998a, ApJ 495, 44+ 
Pettini, M., Katz, N. 



1997, preprint, astro-ph 9712074 



Hernquis t, L., 1998b, preprint, astro-ph 9809401 



Croft, R. A. C, Weinberg, D. H., 
Croft, R. A. C, Weinberg, D. H., 
Dodelson, S., Hui, L., & Jaffe, A. 
Fan, X., 1999, AJ 117, 2528 
Feldman, H. A., Kaiser, N., & Peacock, J. A., 1994, ApJ 426, 23 

Fukugita, M., Ichikawa, T., Gunn, J. E., Doi, M., Shimasaku, K., Schneider, D. P., 1996, AJ 111, 1748 
Gaztanaga, E. & Croft, R. A. C, 1998, MNRAS 309, 885 
Gnedin, N. Y., 1998, MNRAS 299, 392 
Gnedin, N. Y. & Hui, L., 1998, MNRAS 296, 44 
Hamilton, A. J. S., 1997a, MNRAS 289, 285 



21 



Hamilton, A. J. S., 1997b, MNRA S 289, 295 



Hamilton, A. J. S., 1999, preprint astro-ph 9905191 



Hernquist, L., Katz, N., Weinberg, D. H., & Miralda-Escude, J., 1996, ApJ 457, L51 
Home, K., 1986, PASP 98, 609 

Hu, E. M., Kim, T.-S., Cowie, L. L., So ngaila, A., & Rauch, M., 1995, AJ 110, 1526+ 

Hui, L., 1 998, Proceedings of the MPA/ESO Conference Evolution of Large Scale Structure, Garching, Ed. Banday and 

Sheth, |astro-ph 9812293 
Hui, L., 1999, ApJ 516, 519 
2000, in preparation 
& Buries, S., 2000, in preparation 
& Gaztanaga, E., 1999, ApJ 519, 622 
& Gnedin, N. Y., 1997, MNRAS 292, 27 
Gncdin, N. Y., & Zhang, Y., 1997, ApJ 486, 599+ 
Kim, T.-S., Hu, E. M., Cowie, L. L., & Songaila, A., 1997, AJ 
Kirkman, D. & Tytler, D., 1997, ApJ 484, 672+ 
Landy, S. D. & Szalay, A. S., 1993, ApJ 412, 64 

Lu, L., Sargent, W. L. W., Womble, D. S., & Takada-Hidai, M., 1996, ApJ 472, 509+ 

McDonald, P., Miralda-Escud/'e , J., Rauch, M., Sargent, W. L. W., Barlow, T. A., Cen, R., & Ostriker, J. P., 1999, ApJ, 



Hui 
Hui 
Hui 
Hui 
Hui 



3 ., 1992b, Numerical Recipes 

Weinberg, D. H., Hernquist, L., Katz, N., Cen, R., & 



submitted, astro-ph 9911196 
Meiksin, A. & White, M., 1999, MNRAS, submitted, |astro-ph 9812129 
Miralda-Escude, J., Cen, R., Ostriker, J. i. P., & Rauch, M., 1996, ApJ 471, 582- 
Muecket, J. P., Petitjean, P., Kates, R. E., & Riediger, R 1996, AfcA 308 , 17 
Nusser, A. & Haehnelt, M., 1999, submitted to MNRAS, |astro-ph 9906406| 
Peebles, P. J. E., 1980, The large-scale structure of the universe, Princeton University Press 
Peterson, B. M., 1997, An introduction to active galactic nuclei, Cambridge University Press 

Phillips, J., Weinberg, D. H., Croft, R. A. C, Katz, N., Hernquist, L., & Pettini, M., 2000, preprint, |astro-ph 0001089 
Press, W. H., Rybicki, G. B., & Hewitt, J. N., 1992a, ApJ 385, 404 
Press, W. H., Rybicki, G. B., & Schneider, D. P., 1993, ApJ 414, 64 
Press, W. H., Teukolsky, S. A., Vetterling, W. T., & Flannery, B. P., 
Rauch, M., 1998, Ann.Rev.Astr.Ap., in press 
Rauch, M., Miralda-Escude, J., Sargent, W. L. W., Barlow, T. A 

Ostriker, J. P., 1997, ApJ 489, 7+ 
Reisenegger, A. & Miralda-Escude, J., 1995, ApJ 449, 476+ 
Rybicki, G. B. & Press, W. H., 1992, ApJ 398, 169 

Scoccimarro, R., Zaldarriaga, M., & Hui, L., 1999, ApJ, submitted, astro-ph 9901099 
Seljak, U., 1998, ApJ 503, 492 
Szapudi, S. & Szalay, A. S., 1998, ApJ 494, L41 

Tegmark, M., Hamilton, A. J. S., Strauss, M. A., Vogeley, M. S., & Szalay, A. S., 1998, ApJ 499, 555 
Tegmark, M., Taylor, A. N, & Heavens, A. F., 1997, ApJ 480, 22+ 
Vogeley, M. S. & Szalay, A. S., 1996, ApJ 465, 34+ 
Vogt, S. S. et al., 1994, Proc. SPIE 2198, 362 

Weinberg, D. H., Croft, R. A. C, , Hernquist, L., Katz, N., & Pettini, M., 1998, preprint, |astro-ph 9810011 
York, D. et al., 1999, Internal SDSS Report 
Zaldarriaga, M., Seljak, U., & Hui, L., 2000, in preparation 
Zhang, Y., Anninos, P., & Norman, M. L., 1995, ApJ 453, L57 

Zheng, W., Kriss, G. A., Telfer, R. C, Grimes, J. P., & Davidsen, A. F., 1998, ApJ 492, 855 
Zuo, L. & Bond, J. R., 1994, ApJ 423, 73 




Fig. 1. — The lower panel shows the transmission e~ T as a function of wavelength u takep from a SCDM simulation, with no noise added. 
The upper panel shows the corresponding (normalized) transmission power spectrum (eq. P]). All subsequent simulated spectra in this paper 
are based on this one, with various levels of noise, contamination, etc added. 
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Fig. 2. — Bottom two panels show the reduced quasar spectrum (signal) with its noise array (eq. pj]), which resembles a high resolution 
(8 km/s FWHM) echelle spectrum with 12 orders, each 50A long^ It is assumed a relative calibration between_thc different orders has 
been attempted. The spectrum is generated from the one in Fig. hi, with suitable noise added as described in § p.2| . Note that the overall 
normalization of the signal or noise is arbitrary but the ratio of the two is not. The dashed line of the bottom panel shows the input continuum 
leve l The second panel from the top shows the recovered mean-transmission for three different choices of the mean-transmission-basis (p; eq. 
Jig): dotted line for a flat model-continuum/mean (p° = 1 only); short-dashed line for a basis of polynomials up to the third order (p° ... 

long-dashed line for a basis that is the same as that for the short-dashed line except that each echelle order is fi tted separately. The solid 
line give s the true mean. The top-panel shows the result of applying the power spectrum estimator P2 (eq. 12a |) with uniform weighting 
(eq. |14|). The symbol AP2/P stands for (P2 — P)/P where P is the true power spectrum. The different lines correspond to different choices 
of the mean-basis, labeled as in the second panel. 
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Fig. 3. — Similar to Fig. y except that a relative calibration between different echelle orders have not been done, and the input continuum is 
taken from continuum-fits to an observed quasar spectrum. The dotted line in the top two panels corresponds to the case where the continuum 
is modeled as flat for each order. The dashed line is where polynomials up to the third order are used to fit for the mean transmission in each 
order. The solid line in the second panel from the top shows the true mean transmission. 
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Fig. 4. — The bottom 2 panels are the same as in Fig. The solid line in the top second panel shows the true continuum level (as opposed 
to the true transmission mean level as before). The dashed line is a continuum-fit to the simulated spectrum. In the top panel, AP/P 

denotes (P U n — Pun) /Pun (the un-normalized power spectrum eq. [ |lo| fe [EJ) for thc_upper long-dashed line; while it denotes (Pi — P)/P (the 
normalized power spectrum, but estimated using the continuum-fitted data; eq. [nil fe pf) for the lower short-dashed line. 
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Fig. 5. — Similar to Fig. |i| except that the mean transmission is lower (0.39 instead of 0.64). The solid and dashed line in the second 
panel from the top represent respectively the true continuum level and the continuum-fits. In the top panel, the upper long-dashed line 
shows (Pun — Pun)/Pun (error for the un-normalized power spectrum estimated from eq, ^o|); the lower short-dashed line shows (Pi — P)/P 
(error for the normalized power spectrum estimated using continuum-fitted data; eq. pi]); tine dotted line shows (P2 — P)/P (error for the 
normalized power estimated using trend-removal with a fiat trend for each echelle order; eq. JL2|). 
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0FlG. 6. — The bottom 2 panels show a low resolution (FWHM = 1.17 A), somewhat noisy simulated spectrum based on the one in Fig. 
The solid line in the second panel from the top is the true continuum (not mean transmission) whereas the dashed line corresponds to a 
continuum-fit to the data in the bottom panel. In the top panel, AP/P denotes (P U n — Pun) /Pun (the un-normalized power spectrum eq. Jlp[ 
& Jlj]) for the lowpr long-dashed line; while it denotes (Pi — P)/P (the normalized power spectrum, but estimated using the continuum-fitted 
data; eq. [fll| & fcj]) for the upper short-dashed line. 
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Fig. 7. — The bottom 2 panels are the same as in Fig. M. The second panel from the top shows the recovery of the mean transmission 
assuming a model composing of polynomials up to the thirdorder (dashed line). The solid line shows the true mean. The dashed line in the 
top panel represents (P% — P)/P — the error for the normalized power spectrum estimated using eq. (fl^). The dotted line shows the power 
spectrum estimate if shot-noise were not subtracted. 
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Fig. 8. — The bottom two panels show continuum-fits taken from the red (bottom) and blue (second from bottom) sicics^of Lya emission 
from the spectrum of a quasar at redshift 2 = 3. The top panel shows the power spectra of the red and blue continua (eq. 1 21b |) with solid and 
dotted lines respectively. They are very close to each other the panel below it shows their fractional difference [Pc(blue) — Fq, (red)]/Pc (red). 
Also shown in the top panel is the quantity ee7(fc) (eq. |E2|), with solid and dotted lines denoting its values on the red and blue sides. They 
differ by only a few percent. 
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Fig. 9. — A demonstration of power spectrum corrections on large scales. The theoretical spectrum from Fig. gj is multiplied by a set of 
10 different continua (one ofjvhich is shown in the middle panel), convolved with a Gaussian of 0.125 AFWHM and with a small amount of 
noise added (similar to Fig. Eh. The resulting 10 simulated spectra (one of which is shown in the bottom panel) is analyzed and the resulting 
power spectrum fractional error is shown in the top panel. Solid line shows error in the power spectrum estimate with no corrections applied 
(cq. 1 12| or PQI ); dotted line shows the error using the estimator P4 in eq. (p4|), and dashed line shows the error using P3 from eq. (p3|). 
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Fig. 10. — Similar to Fig. H, except that the noise level and resolution resemble instead those of Fig. |6| Note how the noise added makes 
the continuum (middle panel), from which we estimate the continuum power spectrum, quitfi_moisv_as well. In the top panel, the solid line 
shows the fractional error in the pow er spectrum estimate with no corrections applied (eq. [113] or |20| ) , and the dotted line shows the error 
using the estimator P4 from eq. (124). 
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Fig. 11. — The bottom panel shows a simulated spectrum with resolution and S/N similar to Fig. H, except that metal absorption lines 
shown in the second panel have been added. The y-axis of the second panel is e~ T where r is the opticaTdepth due to metal absorption. The 
seond panel on the top shows the true mean transmission and the recovered mean transmission assuming a flat trend for each order (dotted 
line). The dotted line in the top panel is (P2 — P)/P for the case where no attempt is made to cut out the metal lines; the dashed (solid) line 
is the same fractional error for the normalized power spectrum for the case where all metal lines with r > 1 (r > 0.4) are discarded and the 
corresponding gaps filled in via interpolation. 
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Fig. 12. — Power spectrum measurements from low resolution and noisy quasar spectra (similar to Fig. Q), with (dotted) and without 
(dashed) metal contamination. 
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IG. 13. — Power spectrum estimation using uniform weighting (bottom panel; eq. ]l4) ) versus minimum variance weighting (top panel; eq. 
1). The simulated QSO spectrum consists of 12 segments, half of which have comparable S/N to Fig. H and half of which have ~ 20 times 
lower S/N. 
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Fig. 14. — The solid line is the mean observed transmission power spectrum at z = 3 from McDonald ct al. (1999). Its continuation as a 
dashed line at small k's indicate the scales at which continuum might introduce significant uncertainties, and so the power is not shown. The two 
dotted lines at top and bottom shows the level of shot-noise at two extremes. The top one resembles the quality expected for SDSS spectra while 
the bottom resembles HIRES Keck spectra of bright quasars. The interval of Afc uac at the bottom indicates the range of scales that arc currently 
used to recover the linear mass power spectrum. The long tickmarks at the top indicates the resolution of SDSS and typical Keck HIRES spec- 
tra. Different units of k are related by [k/ (km/s)- 1 ] = [(l+z)/4][fc/A _1 ]/61.67 = [(l+2)/y / f2 m (l + z) 3 + Q k (l + z) 2 + nX][fc/hMpc _1 ]/100. 



